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SECTION  1 
INTRODUCTION 


The  University  of  Dayton  Research  Institute  C.UDRI )  has 
completed  an  exploratory  development  contract  with  the  Plight 
Dynamics  Laboratory  (AFWAL/FIBR)  to  develop  a  computer  program 
for  the  design  and  analysis  of  damage  resistant  structures.  In 
addition,  the  Vought  Corporation  utilized  UDRI  computer  output 
to  redesign  an  A7D  aircraft  outer  wing  panel.  This  report 
documents  the  theoretical  development  behind  the  computer  code 
referred  to  as  ADDRESS  {Automated  Design  of  Damage  Resistant 
Structures)  and  presents  applications  to  typical  aircraft 
structures.  Volume  II  of  the  report  gives  details  of  the  computer 
code.  Detailed  design  information  and  a  test  plan  for  the  A7D 
panel  are  submitted  under  separate  cover  as  required  by  the 
contract. 

This  introductory  section  gives  the  background  of  the 
technological  requirement,  and  a  statement  of  the  overall  objective. 
It  concludes  with  a  summary  of  the  major  accomplishments. 

1.1  BACKGROUND 

By  the  very  nature  of  their  reason  for  existence,  military 

aircraft  are  susceptible  to  damage.  The  performance  of  a  damaged 

aircraft  can  be  affected  in  a  variety  of  ways  depending  on  which 

of  the  aircraft  systems  are  affected  and  to  what  degree  the  damage 

is  inflicted.  Structural  damage  to  aircraft  has  been  responsible 
for  a  significant  number  of  aircraft  losses.  Consequently,  it  is 
extremely  important  that  military  aircraft  maintain  an  adequate 
degree  of  structural  integrity  under  damaged  conditions.  The  best 
way  to  insure  the  structural  adequacy  of  damaged  aircraft  is  to 
make  the  proper  provisions  during  the  design  phase. 

One  obvious  way  to  provide  adequate  residual  strength  for 
damaged  structures  is  to  increase  the  size  of  the  various  members 
and  components  which  make  up  the  structures.  This  simplistic 
approach,  however,  is  at  odds  with  the  constant  quest  to  reduce 
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aircraft  weight  and  cost.  Consequently,  aerospace  engineers  strive 
to  design  aircraft  which  are  structurally  adequate  but  which  have 
lower  weight  and/or  cost  than  alternative  designs.  Conventional 
methods  (largely  intuitive)  permit  engineers  to  reduce  the 
weight  of  initial  designs  by  resizing.  However,  these  procedures 
do  not,  in  general,  yield  designs  which  are  "best"  in  any  sense; 
they  are  just  better  than  previous  trials.  In  addition,  vulner¬ 
ability  is  usually  considered  only  as  an  afterthought.  Con¬ 
sequently,  there  is  a  real  need  for  the  development  of  logical 
procedures,  based  on  sound  mathematical  principles,  which  will 
produce  damage  tolerant  structural  designs  that  satisfy  a 
criterion  of  merit  such  as  minimum  weight. 

Although  there  has  been  considerable  progress  in  structural 
optimization,  these  developments  have  been  largely  academic. 

However,  there  is  now  sufficient  evidence  to  indicate  that  optimiza¬ 
tion  can  be  a  valuable  practical  tool  in  reducing  weight  and  cost, 
and  in  improving  performance.  Therefore,  in  view  of  the  above 
mentioned  need  to  provide  damage  tolerant  military  structures,  the 
aims  of  this  current  program  have  been  to  put  optimization  into  the 
damage  tolerant  design  cycle.  The  objective  of  the  program  is 
stated  in  the  following  paragraph. 

1.2  OBJECTIVE 

The  primary  objective  of  this  research  program  was  to 
transition  some  of  the  theoretical  developments  into  structural 
optimization  of  the  last  20  years  into  practical  design  tools  so 
that  significant  cost  reductions  can  be  achieved  in  the  design 
of  future  aerospace  vehicles.  Particular  emphasis  was  focused 
on  damage- susceptible  structures. 

1 . 3  SUMMARY 

An  optimality  criterion  approach  is  presented  which  will  aid 
in  the  design  of  damage  tolerant  structures  subject  to  stress, 
deflection,  and  frequency  requirements.  Damage  conditions  are 
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treated  in  an  integral  manner  in  the  resizing  algorithm.  An 
iterative  reanalysis  procedure  is  used  to  improve  the  efficiency 
of  the  static  analyses  that  are  needed  as  the  optimization 
proceeds.  Several  representative  structural  components  have  been 
resized  with  the  procedure  and  a  test  plan  is  outlined  which  will 
help  to  verify  the  damage  tolerant  characteristics.  For  wing 
structures,  size  and  location  of  damage  are  the  driving  features 
of  the  process.  Criteria  are  presented  to  aid  the  analyst  in 
determining  the  critical  stress  allowables. 

The  major  accomplishments  of  this  contract  are:  { i )  an 
extensive  revision  and  consolidation  of  several  important  Air 
Force  structural  analysis  and  optimization  codes;  (2)  an  in-depth 
investigation  of  structural  reanalysis  techniques;  and  (3)  the 
transfer  of  technology  from  a  research  environment  to  one  of 
design,  fabrication,  and  testing.  Additionally,  the  study  of 
representative  structures  aided  in  the  establishment  of  ballistic 
damage  criteria. 

In  the  computer  code  development  area,  it  is  now  possible  to 
analyze,  resize,  and  assess  the  effect  of  ballistic  damage  on  wing 
structures.  Stresses,  deflections,  and  vibration  characteristics 
can  be  computed  and  damage  sensitivities  established  for  both  metal 
and  composite  structures  within  the  framework  of  a  single  computer 
code.  The  program  has  detailed  documentation  and  can  be  extended 
to  consider  additional  design  requirements. 

Reanalysis  techniques  were  studied  and  mathematical  criteria 
established  for  the  convergence  of  the  static  iterations.  A  new 
design  scaling  technique  was  developed  to  meet  a  fixed  frequency 
requirement  and  an  energy  criterion  simultaneously.  The  efficiency 
of  the  reanalysis  was  found  to  be  highly  dependent  on  the  program¬ 
ming  techniques  used. 

The  redesign  of  the  A7D  outer  wing  panel  presented  herein 
required  that  structural  optimization  technology  be  transferred 
to  an  aircraft  design  environment.  As  a  result  of  the  computer 
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optimization  studies,  aircraft  designers  were  able  to  develop  a 
lighter  weight  and  more  damage  tolerant  wing  structure  than  could 
be  done  by  conventional  means.  Full-scale  detailed  design 
drawings  were  made  and  a  test  plan  developed  to  verify  the 
results  of  the  structural  optimization  code. 


/ 
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SECTION  2 

REVIEW  OF  THE  LITERATURE 

Structural  optimization  algorithms  have  been  continuously 
improved  over  the  last  20  years.  There  are  a  number  of  practical 
approaches  now  available  to  large  scale  structural  optimization 
of  airframe  structures.1  However,  there  has  been  little  attempt 
to  verify  the  validity  of  the  optimized  designs  under  actual 
service  conditions.  Of  particular  concern  for  military  aircraft 
is  the  sensitivity  of  the  optimized  structure  to  loss  or  reduction 
of  stiffness  in  members  due  to  damage.  In  this  section 
optimization  for  damage  tolerance  is  reviewed,  followed  by  a 
discussion  of  the  literature  on  structural  reanalysis. 


2.1  OPTIMIZATION  FOR  DAMAGE  TOLERANCE 

The  body  of  the  literature  on  structural  optimization  is 

now  quite  large,  but  there  have  been  relatively  few  investigations 

dealing  specifically  with  the  optimal  design  of  damage  tolerant 

2 

structures.  Arora,  et  al.  have  included  strength,  stiffness,  and 
damage  requirements  in  an  optimal  design  strategy.  A  design 
sensitivity  analysis  is  coupled  with  a  steepest  descent  algorithm 
to  design  damage  tolerant  structures.  An  application  is  made  to 
an  open  truss  helicopter  tail-boom.  Venkayya,  et  al,1  present 
an  iterative  technique  which  demonstrates  the  feasibility  of 
including  damage  considerations  in  a  design  procedure  based  on  an 
optimality  criterion  approach.  Applications  are  made  to  truss 
and  wing  structures. 

Reference  2  enumerates  an  important  list  of  basic 
definitions  regarding  damage  tolerant  design.  These  definitions 
are  reproduced  here  because  they  apply  to  the  current  work  also. 
They  are  as  follows: 

CD  Damage-tolerant  structure.  A  structure  is  called  damage 
tolerant  or  fail-safe  if  it  continues  to  perform  its  basic 
function  even  after  it  sustains  a  specified  level  of  damage. 
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(2)  Damage  condition.  A  damage  condition  for  the  structure 
is  defined  to  consist  of  complete  or  partial  removal  of  a  selected 
member  or  groups  of  members  of  the  structure.  Some  nodes  of  the 
structure  may  be  removed  as  a  result  of  the  damage. 

(3)  Optimal  damage- tolerant  structure .  A  damage-tolerant 

structure  is  called  optimal  if  it  is  designed  to  minimize  a  merit 
function,  subject  to  constraints  that  must  hold  for  the  undamaged  / 

structure  and  under  projected  damage  conditions. 

The  open  truss  helicopter  tail  boom  problem  studied  in 
Reference  2  contains  108  axial  bar  elements  and  38  nodes.  The 
members  are  grouped  so  that  only  42  are  considered  free  design 
variables.  It  is  important  to  note  that  use  was  not  made  in  this 
study  of  efficient  reanalysis  techniques  to  compute  the  response 
for  the  damage  cases. 

Reference  3  makes  extensive  use  of  reanalysis  procedures 
to  compute  response  to  structural  damage.  Elastic  strength, 
yielding,  and  frequency  and  mode  shape  changes  are  estimated. 

Applications  are  to  optimized  truss  and  wing  structures;  no 
attempt  was  made  to  include  the  reanalysis  procedures  within  the 
framework  of  an  optimization  scheme.  It  was  concluded  that  the 
reanalysis  procedures  were  simple,  accurate,  and  a  natural 
extension  of  the  finite-element  based  structural  optimization 
procedures . 

2.2  REANALYSIS  METHODOLOGY 

In  order  to  perform  an  optimization  of  a  structural  design 
which  is  subject  to  multiple  damage  conditions,  it  is  essential  to 
perform  efficiently  the  reanalysis  of  the  damage  conditions.  This 
subsection  contains  a  review  of  the  literature  which  has  application 
to  reanalysis  problems.  First,  the  methodologies  for  static 
reanalysis  are  considered,  followed  by  techniques  for  predicting 
mode  shape  and  frequency  trends  as  structural  modifications  are 
made. 
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2.2.1  Reanalysis  of  Static  Cases 


Several  categories  of  methods  are  available  for  the 

reanalysis  of  structures  with  local  modifications.  The  methods 

all  involve  the  solution  of  the  modified  linear  equations  and 

may  be  classified  as  follows:  (.1)  direct4  7,  (2)  iterative^'8  , 

(3)  implicit  differentiation"*^  {.4)  reduced  basis^"4  or 

1 7“  1 8 

(5)  mixed  .  In  the  current  study,  methods  one  and  two  were 

investigated  in  detail  due  to  their  relative  simplicity  and 
ease  of  programming. 

By  direct  reanalysis  methods  the  solution  to  the 
modified  static  equations  are  obtained  in  one  step.  References  4, 
5,  and  7  solve  the  modified  equations 

(IK]  -  IK1 ] ) {x>  =  {F}  (2.1) 

based  on  the  stiffness  matrix  of  the  unmodified  system 

[KHxJ  =  If)  (2.2) 


where  [K]  is  the  system  stiffness  matrix,  [K^]  is  the  modification, 

{ F 1  is  the  applied  load,  and  {x}  and  {x}  are  the  response  vectors 

for  the  unmodified  and  modified  structures,  respectively.  Use  of 

5  ... 

the  Sherman-Morrison  identity  ,  a  technique  which  gives  an 
explicit  solution  to  equation  (2.1)  when  [K^J  is  a  vector  outer 
product,  yields 


lx)  = 


IK] 


-1 


R 

-  I 

3=1 


[K]"1{ujHuj}T 


IK] 


-1 


{F } 


(2.3) 


where 


Tj  =  “  °j  1)  1  and 

Oj  and  Uj  are  the  eigenvalues  and  eigenvectors  of  the  modification 

matrix  [K^] ,  and  R  is  the  rank  of  [K^] , 

7 

Kavlie  and  Powell  point  out  the  fact  that  the 
Sherman-Morrison  identify  seems  to  require  the  explicit  inverse 
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of  the  [K]  matrix.  However,  when  equation  (2.2)  is  used  in  (2.3), 
it  can  be  shown  that  simplifications  occur  which  negate  the 
requirement  of  an  explicit  inverse  of  the  stiffness  matrix.  A 
main  drawback  of  the  method  is  more  obviously  the  requirement  of 
the  outer  product  representation  of  the  modification  matrix  in 
terms  of  all  of  its  eigenvectors.  The  eigenanalysis  plus  the 
matrix  multiplications  of  equation  (2.3)  could  easily  exceed  the 
cost  of  resolution  of  equation  (2.1)  by  the  Choleski  method,  for 
example. 

Reanalysis  methods  which  use  iteration  may  be  slow 

to  converge,  if  at  all,  when  large  changes  in  the  design  are  made. 

8 

Kirsch  and  Rubinstein  develop  a  method  to  improve  convergence 
based  on  an  expression  of  the  matrix  of  changes  in  stiffness  as  a 
linear  combination  of  two  matrices.  Phansalker^  makes  similar 
conclusions,  as  does  Switzky^®.  Reference  10  improves  the  conver¬ 
gence  of  the  iterations  by  selecting  good  initial  guesses  for  the 

11 

modified  responses.  Venkayya  develops  an  iterative  procedure 
based  on  a  Taylor  series  expansion  of  the  response.  Details  of 
this  procedure  are  given  in  Section  3.3  herein. 

References  12  and  13  propose  using  implicit  differentia¬ 
tion  of  the  governing  equations  to  obtain  response  gradient  infor¬ 
mation.  In  Reference  14,  a  combined  Taylor  series-iterative 
technique  using  first-order  sensitivities  is  developed  and  applied 
with  good  results  to  space  truss  structures.  Reference  15  documents 
a  similar  study,  and  shows  how  sensitivity  information  can  be 
incorporated  into  an  optimality  criterion  or  a  mathematical 
programming  method  for  structural  optimization. 

1 7-18 

The  work  of  Noor  and  Lowder  develops  a  "mixed 

method"  for  reanalysis;  i.e. ,  the  fundamental  unknowns  include  both 
force  and  displacement  parameters.  The  merits  of  the  Taylor  series 
approach  are  assessed  as  applied  with  the  mixed  approach  and  compared 
with  the  displacement  method.  It  is  also  shown  that  the  displace¬ 
ment  approach  can  be  greatly  improved  by  choosing  the  design 
variables  to  be  reciprocals  of  the  sizing  variables. 
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8 

The  statement  by  Kirsch  and  Rubinstein  that  no  single 
method  of  reanalysis  is  superior  in  all  cases  for  static  problems 
seems  to  summarize  the  reanalysis  state  of  the  art.  A  method 
should  be  evaluated  and  selected  on  the  basis  of  its  ease  of  use, 
computational  efficiency,  and  convergence  characteristics.  In  the 
work  reported  herein,  simple  iterations  in  conjunction  with  an 
acceleration  technique  fit  the  above  criteria. 

2.2.2  Reanalysis  of  Eigenvalue  Problems 

Reanalysis  techniques  for  eigenvalue  problems  are 

.  .  .  .  19-21 

classified  herein  as  either  (1)  derivative-based  or  (2) 

.  22-24 

iterative-based  .  In  the  former  case,  implicit  differentiation 

of  the  eigenvalue  problem  is  performed  to  obtain  first  order  gradient 
information.  The  iterative  techniques  may  in  some  cases  also  lead 
to  gradients.  However,  iterations  can  also  be  used  to  estimate 
discrete  changes  in  the  eigenvalues  and  eigenvectors. 

19 

Fox  and  Kapoor  developed  expressions  for  the  rate  of 

change  of  eigenvalues  and  eigenvectors  of  a  discrete  model  of  a 

vibrating  system.  For  a  given  increment  in  the  structural 

variables,  the  modified  eigenvalues  and  eigenvectors  were  computed 

using  a  one  term  Taylor  series  expansion.  The  rates  of  change  of 
'tlr 

the  i —  eigenvalues  were  expressed  by  Fox  and  Kapoor  in  terms 

of  the  rates  of  change  of  the  mass  and  stiffness  matrices  and  the 
.  th 

i  eigenvector. 

Reference  19  presents  two  formulations  for  the  rates  of 

change  of  the  eigenvectors.  In  the  first  formulation,  the 

t  h 

expression  for  the  i^  generalized  mass  is  differentiated  with 
respect  to  the  variables  and  the  resulting  system  is  solved  for  the 
rate  of  change  of  the  eigenvector.  The  system  of  equations 
involves  only  the  ifc^  vector.  The  second  formulation  uses  an 
eigenvector  expansion  theorem  to  represent  the  eigenvector 
derivative  as  a  linear  combination  of  all  the  vectors  of  the  given 
problem.  The  coefficients  are  then  determined  using  the 
orthogonality  properties  of  the  vectors.  The  disadvantage  of  the 
second  method  is  that  it  requires  that  all  the  vectors  be  computed. 
Fox  and  Kapoor  briefly  comment  that  it  may  be  possible  to  use  a 
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partial  sum  to  represent  the  eigenvector  derivative.  Hence,  for 
small  structural  changes,  not  all  the  vectors  may  be  needed  to 
approximate  the  eigenvector  derivatives. 

Rogers20  explored  the  work  of  Fox  and  Kapoor  further 
and  concluded  that  it  could  also  be  extended  to  eigenvalue 
problems  where  the  matrices  are  non symmetric  and  complex.  The 
expressions  for  the  eigenvalues  and  eigenvectors  included  products 
of  both  the  left  and  right  eigenvectors.  Thus,  it  is  necessary  to 
solve  an  additional  system  of  equations  to  obtain  the  reguired 
derivatives.  The  Rogers  expressions  are  useful  in  studying  design 
trends  in  nonconservative  systems. 

2 1 

Nelson  developed  a  simplified  procedure  for 
determining  eigenvector  derivatives  which  is  more  efficient  than 
previous  methods  when  large  systems  are  considered.  The  matrix 
technique  for  determining  the  eigenvector  derivative  is  similar  to 
the  first  formulation  of  Fox  and  Kapoor  discussed  above.  However, 
the  original  eigensystem  of  rank  (n-1)  is  modified  to  convert  it 
to  rank  n  and  then  solved  for  a  vector  which  is  used  to  compute 
the  eigenvector  derivative.  The  system  of  equations  was  solved 
without  extensive  matrix  manipulations  as  indicated  by  Fox  and 
Kapoor . 

Another  class  of  methods22-24  involves  matrix  iterations. 

2  7 

Rudisill  and  Chu  proposed  an  iteration  based  on  a  rearrangement 

of  terms  in  the  Rayleigh  quotient.  Both  self-adjoint  and  nonself- 

2  3 

adjoint  systems  can  be  handled  in  this  manner.  Andrew  proved 
the  convergence  of  the  Rudisill-Chu  iteration  and  examined  its  rate 
of  convergence.  It  was  concluded  by  Andrew  that  the  iterative 
procedures  should  be  used  rather  than  the  "direct  method"  {i.e. , 
the  approach  in  References  19-21)  when  the  number  of  iterations 
is  less  than  n/3,  where  n  is  the  order  of  the  system. 

Hemming,  Venkayya,  and  Eastep24  studied  the  problem  of 
flutter  speed  degradation  of  damaged  optimized  lifting  surfaces. 
Since  a  modal  flutter  analysis  approach  was  used,  it  was  necessary 
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to  determine  the  effect  of  design  perturbations  (i.e.,  damage)  on 
the  natural  frequencies  and  mode  shapes.  The  perturbed  frequencies 
were  obtained  by  the  use  of  the  Rayleigh  quotient.  Mass  and 
stiffnesses  in  this  quotient  were  for  the  damage  structure.  As  a 
first  approximation,  the  modes  for  the  undamaged  case  were  used  to 
compute  the  approximate  triple  products.  The  modal  vectors  were 
improved  by  solving  a  system  of  equations  iteratively,  which  required 
the  decomposition  of  a  different  coefficient  matrix  for  each 
damage  case. 

In  summary,  it  appears  that  iterative  techniques 
are  more  efficient  than  the  Fox  and  Kapoor  type  of  approach. 

However,  in  the  case  when  changes  in  the  structure  are  large, 
the  convergence  of  these  techniques  is  sometimes  questionable. 

The  work  reported  herein  uses  an  iterative  procedure  somewhat 
along  the  same  lines  as  Reference  24  with  the  improvement  that  a 
new  coefficient  matrix  need  not  be  decomposed  for  each  damage 
case .  In  the  cases  of  large  damage ,  it  is  possible  that  stress 
and  deflection  requirements  are  more  important  than  frequency  and 
mode  shape  considerations.  So  the  convergence  of  the  iterations, 
although  an  interesting  mathematical  problem,  may  be  a  moot 
consideration  when  strength  requirements  are  present. 


11 


SECTION  3 

THEORETICAL  DEVELOPMENT 

In  this  section  the  theoretical  development  of  the  analysis 
and  optimization  procedures  for  damage  resistant  structures  are 
given.  The  finite  element  technique  is  discussed,  followed  by 
the  vibration  analysis  algorithm  used  in  this  work.  The  reanalysis 
procedures  are  next  presented,  and  the  section  concludes  with  the 
presentation  of  various  aspects  of  the  design  optimization  strategy 
such  as  design  scaling. 

3.1  STATIC  ANALYSIS 

The  structural  model  for  the  static  analysis  calculations  is 
based  on  the  finite  element  displacement  method. 

The  current  library  includes  the  following  element  types: 

1.  Bar  (axial  force  member) 

2.  Membrane  triangle 

3.  Membrane  quadrilateral 

4.  Shear  panel 

5 .  Beam 

The  local  coordinate  system  for  the  element  types  shown  in  Figure  3.1 
is  the  same  as  in  Reference  25,  The  bar  element  is  a  constant  strain 
line  element  and  is  equivalent  to  a  NASTRAN  rod  element.  Four 
constant  strain  membrane  triangles  are  used  to  form  the  quadrilateral 
element  by  a  merging  process,  followed  by  static  reduction  to 
eliminate  the  interior  node.  The  shear  panel  is  formed  in  a 
manner  similar  to  the  membrane  quadrilateral;  however,  only  shear 
energy  is  considered  in  forming  the  stiffness  matrix.  The  details 
of  the  beam  element  are  included  in  the  next  subsection.  Further 
details  of  the  elements  other  than  the  beam  are  given  in  Reference 
25. 

The  efficiency  of  the  static  analysis  is  highly  dependent  upon 
the  programming  of  the  solution  technique  to  the  equilibrium 
equations 
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(c)  Quadrilateral  or  Shear  Panel  (d)  Beam  Element 


Figure  3.1.  Elements  and  Local  Coordinate  System 
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[K] {x}  =  {F> 


(3-1) 


where  [K]  is  the  global  stiffness  matrix,  {F}  is  the  vector  (or 
vectors)  of  applied  forces,  and  {x}  is  the  solution  vector  (or 
vectors)  of  displacements  and  rotations.  The  Choleski  decomposition 
technique  was  used  for  these  studies. 

The  stiffness  matrix  is  thus  decomposed  into  the  form 


[K]  =  [L]  [L] 


(3.2) 


where  the  components  of  the  lower  triangular  matrix  [L]  are 

i-1  .  ,1/2 


i.  .  =  fa..  -  j  ) 

11  {  11  k=l  lk  > 


(3.3a) 


3i 


i-1 


‘ii  - 1  » 

11  k=l 


1/ 


ik  VK  "u  '  3>i 


(3.3b) 


After  decomposition,  forward 

[L]  (y)  =  (F>  (3.4) 

and  backward  substitutions 

[L]T  {x)  =  {y}  (3.5) 

are  performed  to  obtain  the  solution.  The  efficiency  of  the 

Choleski  method  is  nearly  double  that  of  the  Gaussian  elimination 

25 

procedure  in  0PTSTAT  .  In  both  ADDRESS  and  0PTSTAT  advantage  is 
taken  of  the  special  "skyline"  form  of  the  stiffness  matrix  so 
that  storage  is  minimized. 


3.2  BEAM  BENDING  ELEMENT 

Figure  3.2  shows  the  coordinate  system  and  cross-sectional 
properties  of  the  beam  bending  element.  The  beam  element  is  a 
straight  bar  with  linearly  varing  properties,  and  is  capable  of 
resisting  axial  forces,  bending  moments  about  two  principal  axes, 
and  twisting  moments.  Node  3  is  used  to  determine  the  cross- 
section  orientation  of  the  beam. 
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n  =  Node  Point  (Shear  Center) 
a  =  Mass  Centroid  Location 
b  =  Geometric  Centroid  Location 

Nodal  Degrees  of  Freedom:  u#  v,  w,  0V,  0„,  0„ 

A  X  Z 


Figure  3.2.  Beam  Element  Geometry. 
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The  shear  center,  mass  centroid,  and  geometric  centroid  are 
allowed  to  be  distinct  points.  Any  transverse  force  applied 
through  the  shear  center  causes  no  torsion  of  the  beam.26  When 
the  element  is  twisted,  the  twist  takes  place  around  the  shear 
center  (center  of  twist) .  When  the  cross-section  has  one  axis  of 
symmetry,  the  shear  center  will  be  located  on  that  axis.  For  cases 
when  there  are  two  axes  of  symmetry,  the  shear  center  will 
coincide  with  the  geometric  centroid.  The  unbalance  provided  by 
the  mass  centroid  offset  is  useful  in  modeling  high  aspect  ratio 
lifting  surfaces  for  vibration  problems. 

The  nodal  degrees  of  freedom  (three  displacements  and  three 
rotations)  lead  to  element  matrices  of  order  12,  Terms  in  these 
matrices  are  functions  of  the  following  properties: 

1.  Mass  centroid  location  -  y  ,  z 

^a  a 

2.  Geometric  centroid  location  -  y,  ,  z, 

b  b 

3.  Cross-sectional  area  -  A(x) 

4.  Moments  of  inertia  -  I  (x) ,  I  (x) ,  J(x) 

y  2 

5.  Geometric  unbalance  -  G  (x) ,  G  (x) ,  G  (x) 

y  z  yz 

6.  Mass  moments  of  inertia  -  K  (x) ,  K  (x) 

y  2 

7.  Mass  unbalance  -  S  (x) ,  S  (x) 

y  z 

8.  Elastic  moduli  -  E,  G 

9.  Material  density  -  p (x)  (mass/length)  or 

p*(x,y,z)  (mass/volume) 

The  above  geometrical  and  inertia  properties  are  interrelated  and 
are  defined  by  the  following  integrals: 

r 

j  dA  =  A*x>  (3.6a) 

A 

|  y2  dA  =  Ia(x)  +  y2  A (x)  =  Iz(x)  (3.6b) 

A 
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(3.6c) 


2  2  — 

Z  dA  =  I^fx)  +  Z£  A(x)  =  I^(X) 

A 

y  dA  =  Gy(x)  =  yb  A (x) 

A 

z  dA  =  Gz  (x)  =  A  (x) 

A 

|  yz  dA  =  Gyz(x)  =  yfa  zfa  A(x) 

A 

I  p*  y2  dA  =  p  K  (x)  +  y2  p  A (x)  =  p  K  (x) 
J  z  a  z 

A 

’  2  2  — 
p*  z  dA  =  p  K  (x)  +  Z  p  A  (x)  =  p  K  (x) 

*  y  *  y 

A 

p*  y  dA  =  S  (x)  «  p  y  A(x) 

J  y  a 

A 

p*  z  dA  =  S  (x)  =  p  z  A (x) 

J  z  a 

A 

p*  yz  dA  =  S  (x)  =  p  y  A(x) 
yz  a 

A 


( 3 . 6d) 

(3. 6e) 

( 3 . 6  f ) 

( 3 . 6  q) 

(  3 . 6  h) 

( 3 . 6  i) 

(3.6  j) 
(3.6  k) 


The  x,  y,  z  displacements  based  on  the  Bernoulli  assumptions 
are,  respectively. 


U  =  u  -  yv,x  -  zw,x 
V  =  V  -  Z0 


W  =  w  +  y0^ 


The  strain-displacement  equations  are 


ex  =  u'x  -  *v'xx  -  2W'xx 


(3.7a) 

(3.7b) 

(3.7c) 

(3.8a) 
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Y  -  -z 0 

xy  x,  x 

Y  =  yd 

XZ  2  X,X 


Qse  of  these  equations,  together  with  the  assumptions 
°y  =  °z  =  TyZ  =  0  and  the  stress-strain  relations 


(3.8b) 

(3.8c) 


o  =  e  e 
x  x 

T  =  G  Y 
xy  1  xy 

T  =  G  Y 
yz  1  xz 

leads  to  the  strain  energy  expression 

U  =  i  |  IEA(x)  n,2x  +  EIz(x)  V,J  +  ETy(x)  w ,Jj 


(3.9a) 

(3.9b) 

(3.9c) 


+  2EIyz(x)v'xx  w-xx  +  GJ(x>  x,x 
-  2EA (x)  u,x(ybv,xx  +  2bw,xx) )  dx 
where  GJ(x)  is  the  torsional  stiffness. 


(3.10) 


Likewise ,  the  equations  (3.7)  can  be  used  in  the  kinetic  energy 
expression 


1  (  —2  —2  —2 
T  =  j  p*  +  V  +  W  ]  dV 

V 


to  obtain  for  constant  beam  density, 

L 

T  °  1 


-%  [  (A(x)  (u2+v2+w2)+K  (w)v  2+K  (x)w  2 

*  J  q  *  x  y  x 

+2Kyz(x)v,xW,x  +  (Kz(x)+Ky(x))e2 

mm  m 

-2A(x)u(yav,x  +  zaw,x) 


(3.11) 


+2A(x)0  (y  w  -  z  v) ]  dx 
xua  a 


(3.12) 
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The  displacement  patterns  are  the  following  polynomials  in 
the  variable  5  =  x/L 


u  = 

v  = 


w  = 


0 


X 


(1-5)^  +  Su2 

(1-3C2+2C3)v1 
+  (3?2-2£3)v2 

(1-3?2+2?3}w1 

+  (3C2-2£3)w2 

U-S)exl  +  56 


+ 

+ 

+ 

+ 

x2 


L(?-2^2+^3)vx1 

L(-52H3)vx2 

LU-2£2+£3)wx1 

l(-^2h3)wx2 


(3.13a) 

(3.13b) 

(3.13c) 

(3.13d) 


where  the  subscripts  1  and  2  denote  the  values  at  the  ends  of  the 
beam.  It  should  also  be  noted  that 


v  .  =  0  - 
XI  zx 


w  .  =  -0  . 
xi  yi 


(3.14a) 

(3.14b) 


since  small  rotations  are  assumed.  Using  these  expressions  in 
Equations  (3.10)  and  (3.12)  leads  to  the  quadratic  forms  for  the 
strain  and  kinetic  energies, 


U  =  \  (q}T  [k]  {q} 

t  =  j  {q)T  fq} 

respectively,  where  the  generalized  coordinates  are  ordered 

{q}  =  <u1,v1,w1,exl,eyl,e2l, 

u2'v2'w2'ex2'ey2'0z2 1 

The  [k]  and  [m]  stiffness  and  mass  matrices  are  of  order  12  and 
include  the  effect  of  a  linear  variation  material  and  geometric 
properties . 


(3.14) 

(3.15) 
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3.3  STATIC  REANALYSIS 

The  static,  elastic  deflection  analysis  for  an  unmodified 
structure  involves  computing  the  solution  of  equations  (3.1)  as 
discussed  in  Section  3.1.  Efficient  reanalysis  methods  should  make 
maximum  use  of  the  Choleski  decomposition  of  the  unmodified  stiffness 
matrix . 

When  a  set  of  structural  elements  are  removed  or  partially 
removed  due  to  damage,  the  reanalysis  problem  becomes 

IK]  {x>  =  {F}  (3.16) 

where  the  modified  stiffness  and  loading  matrices  can  be  written  as 
the  original  quantities  minus  known  perturbations  as  follows: 

[K]  =  [K]  -  [K]_]  (3.17) 

{F}  =  {F}  -  {Fx}  (3.18) 


The  modified  stiffness  matrix  can  be  positive  definite  for  small 
damage;  however,  IK^ ]  will  be  just  positive  since  it  contains  at 
least  n-m  rows  and  columns  of  zeros  where  m  is  the  total  number  of 
degrees  of  freedom  affected  by  the  damaged  elements. 

3 

The  work  by  Venkayya,  et  al.  shows  how  the  perturbed  response 
vector,  {x},  can  be  iteratively  determined  without  directly  solving 
equation  (3.16).  The  reanalysis  formula  which  was  developed  is 

{x(V+1)}  =  [A]  {x(V)}  +  {x}  (3.19) 

where  [A]  =  [K J  ^  [K^] .  This  iteration  starts  with  {x^0*}  =  {x}. 

To  derive  an  iteration  scheme  including  {f.^},  introduce  the 
perturbation  parameter,  e,  which  scales  the  damage  terms  of  the 
stiffness  matrix  and  load  vector  as  follows: 


[Kj.]  =  e[Ko) 


(3.20a) 


{  f2  }  =  e{F0} 


(3.20b) 
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As  in  classical  perturbation  theory,  expand  the  modified  response 
in  a  power  series 


(x>  =  {xQ}  +  efx^}  +  £2{x2)  +  e3{x3>  +  ...  (3.21) 

where  {xq}  -  {x}  and  {x^},  {x2),  ...  are  undetermined  vectors. 

Next,  substitute  equations  (3.20)  and  (3.21)  into  equation 
(3.16)  and  group  terms  into  coefficients  of  e.  This  gives  a 
homogeneous  system  which  implies  that  all  coefficients  are  zero. 

The  following  system  is  thus  obtained: 


<Y0>  =  {x>  (3.22a) 

{y-L }  =  [A]{yo>  -  (K]"1{P1}  (3.22b) 

(yi>  =  [A]{yi_1}  i  =  2,3,...  (3.22c) 

where  f y ^ }  -  e^{x^}  and  [A]  =  [K]  A  (v+1) ^-order  accurate 

estimate  of  {x}  is  obtained  by  adding  the  first  v  equations  of 
(3.22)  to  get 


{xV+1}  =  IA]{x{v)}  +  {x}  -  [K]  1  {F]_}  (3.23) 

It  is  important  to  note  that  [K]  ^  is  never  formed  when  equation 
(3.23)  is  applied.  The  next  response  iterate  is  determined  by 
forward  and  back  substitution  based  on  the  decomposition  of  the 
unmodified  stiffness  matrix. 


In  the  case  that  {F^  =  {0},  equation  (3.23)  has  the  form 


{x(V+1)}  =  f[A]V+1  +  [A]V  +  ...  +  [A] + [I ) 


(x)  (3.24) 


It  is  easy  to  show  that  equation  (3.24)  is  the  first  V+1  terms  in 
a  Taylor  series  in  e  as  follows.  Note  that  the  Taylor  series 
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(3.25) 


(x(e) } 


(x)  +  e {x^ (0) } 


2  _  v+1 

+  jy{x2 (0) }+. . -+(v+1j , 1XV+1 (0) } 


+  •  *  * 


where  {x,.  (0)  }  denotes  the  j—  derivative  of  {x}  at  e=0.  By  direct 
differentiation  of  equation  (3.16)  it  can  be  shown  that 

{x.  (0)}  =  n[K]-1{x.  .(0)  >  =  ii-  lA]n{x}  (3.26) 

3  ^  j 

e  J 

Use  of  equation  (3.26)  in  (3.25)  yields  equation  (3.24).  Thus  each 
cycle  of  simple  iteration  is  equivalent  to  adding  an  additional 
higher  order  term  to  the  Taylor  series. 

An  approach  to  a  convergence  study  of  equation  (3.23)  is  most 
easily  demonstrated  by  considering  the  }  =  {0}  case.  The 
equation  (3.24)  is  reminiscent  of  an  iterative  method  for  solving 
a  linear  eigenvalue  problem.  Therefore,  we  can  study  the  conver¬ 
gence  of  equation  (3.24)  by  considering  the  following  eigenvalue 
problem 


y±  =  [A]  {$i}  (3.25) 

th 

where  and  are  the  i  -eigenvalue  and  eigenvector,  respectively, 

or  equivalently, 

Pi(K]{(t>i)  =  tK1Hti}  (3.26) 

Since  [K]  and  [K^J  are  at  least  positive  semi-definite,  with  t K J 
positive  definite,  then  {x}  has  a  unique  expansion  in  terms  of  the 
[K] -orthonormalized  eigenvector  of  [A], 

n 

{x}  =  a  l  c,  {*. }  (3.27) 

i=l  1  1 

fp  m 

where  [K]  (x}/a  and  a  =  max({<jK}  [K]  {x}).  Substitution 

of  equation  (3.27)  into  (3.24)  gives  a  geometric  series  in 
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(v+1) 

lx 


n  v+1 


>  =  l  I  C. 

i=l  j=0  1  1  1 


(3.28) 


Now,  if 


Pi  <  1 


(3.29) 


the  series  converges  to 

n 

lx}  =  a  l  (c./(l-y . ) ) U . }  (3.30) 

i=l  1  11 

where  the  eigenvalues  have  been  ordered  as  Ui  _>  y2  _>  ...  _>  y^. 

Thus  equation  (3.30)  is  the  theoretical  limit  of  the  iteration 
formula  (3.24)  subject  to  the  convergence  criteria  of  equation 
(3.29) . 

A  useful  quantity  in  damage  assessment  can  be  the  ratio  of 
the  strain  energy  in  the  damaged  structure  to  that  of  the  undamaged 
one,  i. e. , 


U/U  =  {x}TIK){x}/{x}T[K] lx)  (3.31) 

By  using  equation  (3.27)  and  (3.3),  this  is  computed  to  be 

u/u  =  I  <c.2/(l-y.))/  l  c.2  (3.32) 

i=l  1  1  i=l  1 

which  by  the  extremal  properties  of  Rayleigh's  quotient  is  bounded 
by 


i/(i-un)  <  U/U  <  1/ ( 1— y ! )  .  (3.33) 

Thus,  if  we  define  a  sequence  of  stiffness  damage  perturbation 

matrices,  [K.  ],  j=l,...,N,  where  N  is  the  number  of  damage  cases, 

and  we  define  y.  J  to  be  the  maximum  eigenvalue  associated  with 
"til 

the  j  eigenproblem  (see  equation  (3.26)),  then  the  potentially 
most  dangerous  damage  case  has 
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P1*  =  max  (Mi  ) 

j 

as  its  largest  eigenvalue.  It  is  easily  verified  that  a  first 
order  estimate  of  is 

Wj*  =  max{l+{UT(Il]-[K.])  {?})/6  (3.34) 

where  ($,{?})  are  the  dominant  eigenpair  of  the  undamaged  stiffness 
matrix?  i.e.,  (8,{£j})  satisfy 

{K]{£}  =  8(5)  (3.35) 

This  system  is  easily  solved  by  matrix  iteration,  and,  hence, 
a  trial  estimate  of  the  most  critical  damage  case  can  be  quickly 
identified  by  Equation  (3.34). 
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3.3.1  Collapse  Analysis  Based  on  Direct  Stiffness  Variation 


form 


If  the  [K^]  modification  matrix  is  of  the  special 


[K1 }  =  [I]  (3.36) 

the  effect  of  direct  stiffness  variations  (changes  in  the  diagonal 
terms  of  the  stiffness  matrix)  can  be  quickly  assessed  by  the 
techniques  of  this  section.  For  changes  of  the  form  above,  the 
equivalent  iterative  expression  for  the  response  of  the  modified 
system  is 


{^hUo1^  tf  {*iJ 


where  (^,{<Jk))  are  the  eigenpairs  of 


(3.37) 


(3.38) 


In  order  for  the  above  iterations  to  converge, 


|e[  <  ^  i  =  l,...,n  (3,39) 

To  see  how  this  iteration  and  convergence  criterion 
can  be  applied,  consider  the  following  simple  example.  Let 


in 

i — 

4 

3  “ 

y 

4 

7 

4 

,  {F>  = 

i 

_  3 

4 

4  _ 

JL. 

and  e  =  2.  The  first  four  iterates  are 


"  .793 

“  .840  " 

~  .884  " 

.  711 

1.069 

9 

.731 

t 

1.281 

9 

.414 

.  586 

1.095 

.258 

1.673 

The  iterates  continue  to  increase  in  vector  norm  and  essentially 
become  unbounded. 
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An  eigenanalysis  using  the  power  method  on  this  matrix 
shows  after  five  iterations  starting  with  a  vector  with  unit 
entries,  gives  =  12.9,  the  largest  eigenvalue.  Using  the 
eigenshift  technique  on  the  system 

(  [K]  -  [I] )  {*}  =  pi  £ 4> >  (3.40) 

where 


gives  =  1.2  as  the  smallest  eigenvalue.  Again  the  power  method 
is  easily  applied.  Thus,  by  equation  (3.39)  we  see  that  direct 
stiffness  modifications  must  be  less  than  1.2  in  magnitude.  When 
they  are  larger  than  this,  it  causes  the  stiffness  matrix  to  no 
longer  be  positive  definite  and  hence  a  kinematic  instability  or 
collapse  mode  is  possible. 

In  the  above  example,  if  e  is  reduced  to  1.1,  simple 
iteration  will  converge  to  three  place  accuracy  in  67  iterations. 
Clearly,  the  need  exists  to  accelerate  the  convergence  of  the 
iterations  if  possible. 

3.3.2  Acceleration  of  Convergence 

Iterations  of  the  form  given  in  equation  (3.23)  are 

an  n-dimensional  generalization  of  the  simple  fixed  point 

27—  2  8 

iteration  problem  : 


xn+l  =  9<xn)  (3-41> 

This  generates  a  sequence  x^  x2,...  which  we  want 

lim  x  =  £ 
n 

n-*» 

Figure  3.3  shows  a  convergent  fixed  point  iteration  and  accelera- 
tion  by  Aitkin* s-A  process.  From  the  figure,  note  that 
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Figure  3.3.  Acceleration  of  Fixed-Point  Iteration  by  Aitken's-A  Process. 


When  r  <1  and 
n 


,  =  (x  -x  ..)/{x  --x  „ ) 

n-1  n  n-1  f  n-1  n-2 

(3.42) 

r  =  (x  .,-x  )/ (x  -x  .) 
n  n+1  n  '  n  n-1 

(3.43) 

r  -r  , 
n  n-1 


n 


<  £ 


(3.44) 


2 

for  some  small  e  >  0,  then  Aitken's-A  process  consists  of 
extrapolating  the  straight  line  segment  joining  points  2  and  3  to 
the  intersection  with  the  line  y=x  {point  4) .  The  rationale  is 
that  if  r  <  1  and  equation  (3.44)  holds,  then  the  g(x)  curve  slope 
is  not  rapidly  changing  and  then  a  straight  line  extrapolation 
will  intersect  the  line  y=x  near  £.  The  coordinate  for  the 
extrapolated  point  is 

2 


x 

n 


n+1 


(Ax  ) 
n 

72 - 

A  X. 


n-1 


(3.45) 


where 

Axk  ■  xk+i  -  xk 

and 

i2xk  ■  ixk+i  ■  4xk 

In  the  ADDRESS  code,  the  vector  of  modified  responses 
for  each  damage  case  are  computed  by  equation  (3.23).  For  two 
consecutive  iterates,  the  error  is  defined  as 


(v+1) 


=  max 
3=1, . . ,L 


r  N 

[l 

Li=l 


<x.(v+1) 

Xi  j 


21  V2 


(3.46) 


where  L  is  the  number  of  load  cases  and  N  is  the  number  of 
unconstrained  degrees  of  freedom.  The  ratio  of  consecutive  errors 
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(V  +  1) 


e(v+l)/e(v) 


(3.47) 


is  tested  to  make  sure  it  is  less  than  1.0  and 


(r(v+i)  _  r(W)/r(v+l)  <  ml 


(3.48) 


is  required  before  the  extrapolation 


-  (v+l) 

x .  .  =  x. 

13  lj 


(x  Cv+1)  -  X  {V) 

ixij  xij 

(v+i)  (vTT 

x.  .  -2x .  .  +x. 


Tv-1) 


13 


13 


13 


(3.49) 


is  performed.  Equation  (3.49)  is  the  vector  analog  of  equation  (3.45). 
Equations  (3.44)  and  (3.48)  are  analogous  where  e  is  selected  to  be  0.1. 

3.4  VIBRATION  ANALYSIS 

This  subsection  describes  the  analytical  technique  that  was 
used  to  solve  the  eigenvalue  problem 

* [M] {x}  =  IK] {x}  (3.50) 

where  IM]  is  the  global  mass  matrix  (including  non-structural  mass 
items) ,  and  IK]  is  the  global  stiffness  matrix.  The  solution  to 
equation  (3.5  0)  is  given  by  U^Ix^}),  i=l,...,n  where  is  the  itlri 
eigenvalue  (natural  frequency  squared),  and  {x^}  is  the  i^  eigen¬ 
vector  (mode  shape).  The  order  of  the  system  is  n,  and  the  lowest  m 
(m  <  n)  frequencies  and  corresponding  mode  shapes,  and  the  generalized 
masses  are  desired.  The  generalized  masses  are  defined  as 

mi  =  {xi^T  IMHxi>  i  =  1, .  . .  ,m<n  (3.51) 

The  ADDRESS  system  uses  a  simultaneous  iteration  algorithm 
2  9 

(SIMIT)  to  compute  the  modes  and  frequencies.  This  method  was 
selected  because  of  its  excellent  convergence  characteristics  and  its 
relative  simplicity.  The  procedure  takes  a  minimum  of  storage  and  is 
easily  adapted  for  the  calculation  of  unrestrained  modes.  The  basic 
steps  of  SIMIT  are  as  follows: 
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1.  Set  random  trial  vectors  [U]  and  orthogonalize 


2.  Solve  for  [Y] :  lL]T[Y]  =  IUJ 


3.  Form  [Z]s  [Z]  =  IM]  [Y] 

4.  Solve  for  [V] :  [L]  [V]  =  [Z] 

5.  Form  [W] :  [W]  =  [V]  [T*]  and  estimate  A^  (see  below 

for  definition  of  IT*]) 

6.  Form  [U]  by  Schmidt  orthogonalization  of  [W] 

7.  Test  |!u-U]j  <  e 

8.  If  test  of  Step  7  passes,  then  solve  [L]T  [X]  =  [U] 

and  print  the  eigenvectors  IX]  and  eigenvalues 
If  test  of  Step  7  fails,  then  set  [U]  =  lu]  and 
go  to  Step  2 . 

The  key  step  in  the  method  is  a  matrix  iteration  of  the  form 

[A]  [U]  =  [V]  where  [A]  is  a  symmetric  matrix,  [U]  is  an  n  by  m 

matrix  of  iterated  vectors.  Note  that  all  modes  of  interest  are 
simultaneously  iterated.  In  practice,  it  is  necessary  to  use  about 
twice  as  many  trial  vectors  as  the  number  of  vectors  desired.  The 
SIMIT  method  differs  from  the  usual  matrix  iteration  method  in  that 
a  sweeping  matrix  is  not  required.  Also,  the  [ K ] — 1 [M]  matrix  is 
never  formed. 

As  stated  above,  the  matrix  [A]  must  be  symmetric.  An 
important  step  in  the  SIMIT  technique  is  to  transform  the  eigen¬ 
value  problem  of  equation  (3.50)  to  the  form 


[A] {q>  =  A {q } 


(3.52) 


where  [A]  is  symmetric,  which  has  the  eigensolutions  (A  ^ ,  qi> , 
i  =  n.  This  can  be  accomplished  by  defining 


IA]  =  [L]  1  [M]  [L] 
{q>  «  IL]T{x} 


(3.53a) 


(3.53b) 
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where  [K]  =  [L] [L]^.  Note  that  {A]  so  defined  is  symmetric  and 
that  the  eigenvalue  problem  described  by  Equations  {3.52)  and 
{3.53)  is  equivalent  to  equation  (3.50). 

In  practice,  the  matrix  [A]  is  never  explicitly  formed. 
Instead,  for  a  given  set  of  trial  vectors  [U] ,  the  product  I A] [U] 
is  formed  in  the  following  steps: 

1.  Solve  [L]T[Y]  =  [U]  for  [Y]  by  back  substitution 

2.  Form  [2]  =  IM] [Y]  by  premultiplication  of  [Y]  by  [M] 

3.  Solve  IL] [V]  =  [Z]  for  [V]  by  forward  substitution. 

It  is  easily  verified  that  the  above  steps  which  make  use  of  the 
properties  of  the  [Lj  matrix  are  equivalent  to  forming  the  matrix 
product  [a] [U] .  In  the  ADDRESS  program,  IL]  is  stored  in  "skyline" 
or  variable  bandwidth  form. 

The  diagonal  terms  of  [U]T[V]  are  used  as  the  first  estimates 
of  the  eigenvalues.  The  eigenvectors  must  be  modified  before  the 
next  step,  since  they  would  otherwise  all  converge  to  the  vector 
corresponding  to  the  largest  eigenvalue.  In  the  terminology  of 
Reference  29,  an  "interaction  analysis"  is  performed  on  the  vector 
which  decouples  them.  The  decoupled  vectors,  [W] ,  follow  from  the 
matrix  multiplication 


[W]  =  [V] [T* ]  (3.54) 

where  [T*]  has  unit  diagonal  terms  and  skew-symmetric  off-diagonal 
terms  given  by 


t* 


ij 


2b .  ./c .  . 

xy  xj 


(3.55) 


where  b^  are  the  elements  of  the  "interaction  matrix" 

[B]  =  [U]T[V]  (3.56) 

and 
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and 


c 


ij 


b.  .  + 
33 


-b.  . )  2+16b . 

33  iD 


(3.57) 


with  S  equal  to  the  sign  of  b^-b^.  The  [W]  vectors  are  then  made 
orthogonal  by  the  usual  Grahm  Schmidt  orthogonalization  process.  By 
efficient  use  of  storage  registers,  it  is  not  necessary  to  have 
separate  storage  for  the  [B]  and  [T*]  matrices. 

The  orthogonalized  vectors  [U] ,  obtained  from  the  [W]  matrix 
and  then  stored  in  [U]  and  the  iteration  is  continued.  The  process 
is  completed  by  performing  either  a  specific  number  of  iterations 
or  by  achieving  a  given  small  error  on  the  Euclidean  norm  of 

ruj-ru] . 

The  final  step  of  the  procedure  is  to  use  back  substitution 
to  recover  the  mode  shapes  for  the  original  vibration  problem  of 
equation  (3.50).  This  is  accomplished  by  solving  equation  (3.53b) 
where  the  columns  of  the  converged  [U]  are  the  {q^}  vectors. 

A  number  of  investigators  have  used  SIMIT  (see  Reference  29) 
on  large  sized  problems.  Although  the  basis  for  the  interaction 
analysis  as  presented  in  Reference  29  is  somewhat  empirical, 
numerical  results  have  been  consistently  good. 

In  the  case  of  unrestrained  structures,  the  [K]  matrix  will 
not  be  positive  definite.  In  this  case,  use  can  be  made  of  the 
eigenvalue  shifting  technique.  The  quantity  u[M]{x.)  is  added  to 
both  sides  of  equation  (3.50)  where  u  is  a  predetermined  positive 
constant.  Hence,  the  "shifted  problem"  becomes 


[k]{x>  =  nlMHx} 

where 

(3.58a) 

IK]  =  IK]  +  ulM] 

(3.58b) 

n  =  A  +  y 

For  the  work  herein,  jj  is  selected  as 

(3.58c) 
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(3.58d) 


p  =  min  (K . ./M . . ) 
j  33  33 

where  K..  and  M.,  are  the  diagonal  terms  of  the  stiffness  and  mass 
33  33 

matrices,  respectively.  Hence,  for  unrestrained  modal  calculations, 
equation  (3.58a)  is  solved  by  the  SIMIT  technique,  and  the  desired 
eigenvalues  are  then  recovered  by  solving  equation  (3.58c)  for 
It  is  important  to  note  that  even  the  zero  frequency  node  shapes  are 
obtained  in  the  solution  of  equation  (3.58a). 

3.5  FREQUENCY  AND  MODE  SHAPE  REANALYSIS 

A  number  of  papers  have  been  written  on  the  rates  of  change 

of  eigenvalues  and  eigenvectors  with  respect  to  system  parameters 

24 

(see  References  19-23,  for  example) .  Hemming,  et  al.  have 
developed  a  technique  specifically  for  the  analysis  of  damaged 
structures.  The  eigenvalues  are  determined  from  a  Rayleigh  quotient 
procedure  for  the  frequency  calculation;  however,  the  modified 
vectors  are  iteratively  determined  from  solving  a  system  of  equations 
based  on  the  decomposition  of  the  stiffness  matrix  for  the  undamaged 
structure . 

The  reanalysis  procedure  is  derived  by  first  considering  the 
vibration  problem  for  the  undamaged  structure 

[K]{xiJ  =  AilM]{xi}  i  -  1, . . . ,m  (3.59) 

and  a  particular  damage  case 

[KHXi)  =  Xi[M]{5e.}  i  *  1, . . .  ,m  (3.60) 

where 

[M]  -  [M]  -  [M^]  =  modified  mass  matrix 
{x.}  =  {x^}  -  {x,^}  =  modified  eigenvectors 
-  A^  =  modified  eigenvalues 

and  [K]  is  given  by  equation  (3.17). 
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The  iteration  equation  for  the  i —  modified  eigenvalue  has 
been  selected  to  be  the  Rayleigh  quotient  for  equation  (3.60), 

*1iw  -{vv>}T  ^vvy{*i(v)}T®{vv>}  <3-6i) 

where  (v)  denotes  the  iteration  counter.  The  starting  mode  shapes 
on  the  right  side  of  equation  (3.61)  are  taken  to  be  the  mode  shapes 
of  the  undamaged  system.  The  iteration  equations  for  the  ith 
eigenvector  are  obtained  by  substituting  equation  (3.61)  into 
(3.60)  and  rewriting  in  the  iterative  form 

[K]^  (V+1)}  =  (Xi(v)[M]  +  lK1])|xi(V)|  (3.62) 

The  flow  of  calculations  is  then  as  follows: 

1.  Solve  equation  (3.59)  by  simultaneous  iteration  and 
save  the  decomposition  of  the  stiffness  matrix. 

2.  Evaluate  equation  (3.61)  to  get  the  first  estimates 
of  the  eigenvalues • 

3.  Solve  equation  (3.62)  using  the  decomposition  of  the 
stiffness  matrix  determined  in  Step  1.  Iterate  this  equation 
several  times  until  the  norm  of  the  vectors  changes  less  than 
some  specified  amount. 

4.  Orthogonalize  the  last  set  of  iterated  vectors  from 
Step  3  with  respect  to  the  modified  stiffness  matrix. 

5.  Repeat  Steps  2,  3,  and  4  a  prespecified  number  of  times. 
Each  time  Step  2  is  performed,  use  the  last  set  of  eigenvectors 
obtained  from  Step  4  as  input. 

The  iteration  scheme  used  herein  is  similar  to  that  of 

■  .  22 

Rudisill  and  Chu.  In  their  work,  however,  they  are  considering 
the  derivative  problem  rather  than  the  perturbation  problem. 

In  Reference  22,  the  eigenproblem  is  differentiated  and  terms  are 
rearranged  so  that  eigenvector  derivatives  appear  on  both  sides 
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of  the  equation.  This  form  is  then  used  as  the  basis  for  a 
fixed-point  iteration  scheme. 

Several  key  questions  arise  in  the  use  of  iteration 
expressions  such  as  these: 

1.  Do  the  iterations  converge? 

2.  Do  they  converge  to  a  unique  solution? 

3.  Is  it  the  correct  solution? 

Reference  23  deals  with  such  questions  for  the  Rudisill-Chu  itera¬ 
tion  method.  The  method  converges  for  the  dominant  eigenvector/ 
eigenvalue;  however,  in  other  cases  the  solution  characteristics 
can  depend  on  the  initial  guess  for  the  eigenvector  derivative. 

For  the  perturbation  method  used  herein,  satisfactory 
results  have  been  obtained  when  the  mass  and  stiffness  changes 
are  small.  This  implies  that  mode  shape  changes  are  not  large. 
Numerical  results  are  presented  in  Section  5  herein. 

3.6  OPTIMIZATION  FOR  STRESS  AND  DEFLECTION  REQUIREMENTS 

Following  a  development  similar  to  the  optimality  criterion 
approach  of  Reference  30,  the  object  of  the  optimization  process  is 
to  provide  a  structure  which  is  satisfactory  from  a  stress  and 
flection  standpoint  and  has  minimum  total  weight.  It  is  further 
required  that  the  structure  also  satisfy  these  requirements  in  its 
damaged  configuration  or  configurations. 

The  total  structural  weight  of  a  finite  element  model  with 
"m"  elements  is 


m 

W  =--iPi  Ai  li  (3.63) 

where  pi  is  the  mass  density  and  A^  is  volume  of  the  ith  element. 
If  the  elements  are  bars,  then  Ai  is  the  cross-sectional  area  and 
is  the  length.  For  membrane  and  shear  panels  the  A^  and  l . 
are  the  thickness  and  area,  respectively.  The  beam  element 
variables  are  the  same  as  for  the  bars. 
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A  necessary  condition  for  an  optimal  design  is  that  the 
function 

m 

<f>  =  l  P.  A  .A.  -  \Z  {3.64) 

i=l  1  11 

be  stationary  with  respect  to  the  design  variables ,  A^ ,  where  \ 
is  a  Lagrange  multiplier  and  Z  is  a  stiffness  requirement.  When 
Z  is  selected  as  the  strain  energy  of  the  structure,  it  can  be  shown 
that  the  stationary  conditions  lead  to  the  following  statement  of  the 
optimality  criterion  for  the  generalized  stiffness  requirement: 

"The  optimum  structure  for  a  generalized 
stiffness  requirement  is  the  one  in  which  the 
ratio  of  the  average  strain  energy  density  to 
its  mass  density  is  the  same  for  all  elements." 

When  a  displacement  constraint  is  present,  the  stiffness 
requirement  has  the  form  of  the  virtual  work  expression 

Z  =  {F}T  (r>  (3.65) 

where  {F}  is  a  unit  load  vector  and  {r}  is  the  generalized  displace¬ 
ment  vector.  The  unit  load  vector  contains  zeroes  except  in  the 
unit  entry  corresponding  to  the  displacement  constraint.  Again,  from 
the  stationary  condition  of  the  ()>  function,  the  following  statement 
of  the  optimality  criterion  is  given  as  follows: 

"The  optimum  structure  for  a  specified 
displacement  is  the  one  in  which  the  ratio 
of  average  virtual  strain  energy  density  to 
the  mass  density  is  the  same  for  all  elements." 

Damage  cases  are  handled  in  the  same  way  that  multiple  loading 
conditions  are  treated.  The  optimality  conditions  then  take  the 
form 


ed> 

'1  i 


*2ei 


(2) 


+  x  e .  ^  =  l 

pi 


(3.66) 


i  =  1 , . . . ,m 
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where 


e  . 

1 


4  {r}TrK.]  {r.}  /  P.  A.».  (Stress) 

7  l  3/111 

{f .}T[K.) {r.}  /  p.  Aili  (Deflection) 


(3.67) 


th 

and  {f.}  is  the  displacement  vector  corresponding  to  the  j  unit 
load  vector.  The  [K^J  matrices  are  the  stiffness  contributions 
of  the  elements  for  i=l,...,m.  The  index  p  equals  the  number 
of  damage  cases  times  one  plus  the  number  of  loading  conditions. 


Resizing  of  the  element  areas  takes  place  in  accordance 
with  formulas  of  the  type 


(ctiA)v+l  =  (3.68) 

where  cr  is  the  normalized  member  size,  A  is  a  scale  factor,  v 
denotes  the  iteration,  and  E.  is  a  weighted  average  of  the  strain 

4- 1 

energy  density  in  the  1  element.  The  scale  factor  can  be  exactly 
determined  for  stress  and  deflection  constraints. 

In  the  ADDRESS  computer  code,  resizing  takes  place  first  for 

the  stress  requirements  and  then  for  displacements.  A  specified 

number  of  resizings  take  place  to  make  the  strain  energy  density 

uniform  in  each  element  when  subjected  to  the  load  case/damage 

conditions.  The  design  is  scaled  so  that  all  the  stress  and 

deflection  requirements  (if  any)  are  satisfied.  Resizing  is 

terminated  in  the  energy  mode  (i.e.,  stress  constraints)  if  the 

weight  increases.  Resizing  then  takes  place  for  the  displacement 

requirements  and  checks  are  made  on  stress  levels  and  appropriate 

scalings  of  member  sizes  are  made.  The  algorithm  is  very  similar 
25 

to  the  0PTSTAT  program  except  that  additonal  "pseudo  load 
cases"  (i.e,,  damage  cases)  are  included. 
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3.7  OPTIMIZATION  FOR  VIBRATION  REQUIREMENTS 

Reference  31  develops  an  optimality  criterion  approach  to 
the  minimum  weight  design  of  structures  subjected  to  dynamic  loads. 

In  many  design  problems,  it  is  possible  to  characterize  the 
dynamic  load  problem  as  a  requirement  that  a  given  natural  frequency 
not  be  lower  than  some  specified  value.  Such  frequency  constraint 
problems  are  solvable  by  direct  search  techniques  which  involve  the 
calculation  of  frequency  gradients. 

The  optimality  criterion  approach  does  not  treat  a  frequency 
as  a  formal  constraint  as  indicated  above,  but  rather  the  object  is 
to  obtain  minimum  weight  with  maximum  stiffness.  The  optimality 
criterion  is  that  difference  between  the  strain  energy  density  and 
the  kinetic  energy  density  {Lagrangian  density)  be  uniform  for  all 
the  elements  when  the  structure  is  vibrating  in  its  mode  of  interest. 
The  resizing  formula  is  the  case  in  equation  (3.68)  where  E.  is  the 
Lagrangian  density.  The  procedure^  for  using  the  resizing  formula 
is  as  follows: 

1.  With  an  assumed  relative  design  vector  (such  as 

aI  ~~  °2  —  •  ~  =  1.0)  ,  the  desired  normal  mode  and  corresponding 

frequency  are  determined. 

2.  The  relative  strain  energy  and  kinetic  enerqy  of  each 
element  in  this  mode  is  determined. 

3.  Each  design  variable  for  the  next  iteration  is  determined 
by  substituting  the  relative  energies  in  equation  (3.68). 

4.  The  design  vector  is  normalized  with  respect  to  the 
largest  variable.  This  normalization  eliminates  the  need  for 
determining  the  constant  of  proportionality,  A.  (The  next  section  is 
concerned  with  determining  this  constant.) 

5.  The  procedure  is  repeated  until  the  optimality  criterion 
is  attained  or  a  prespecified  number  of  iterations  have  taken  place. 

Experience  with  the  ADDRESS  code  has  indicated  that  in  the 
case  of  damage  with  small  non structural  mass ,  a  frequency  requirement 
is  not  a  meaningful  constraint,  it  is  of  interest,  however,  to 
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check  the  frequencies  before  and  after  damage.  The  frequencies  of 
the  optimal  structure  which  is  designed  for  stress,  deflection,  and 
damage  tolerance  under  multiple  loads  does  not  have  frequencies  which 
are  sensitive  to  changes  in  design  variables. 

If  the  structure  is  designed  from  only  a  frequency  standpoint, 
it  may  be  possible  that  the  particular  mode  being  tracked  is  no 
longer  the  mode  with  the  smallest  eigenvalue.  For  example,  a  uniform 
wing  design  may  have  its  lowest  frequency  in  the  first  bending  mode. 
After  several  resizings,  it  may  be  possible  for  the  first  torsion 
mode  frequency  to  be  below  the  bending  frequency.  At  this  point  the 
algorithm  must  shift  to  computing  energies  for  the  displacement 
pattern  of  the  torsion  mode.  Although  such  conditions  were  not 
noted  in  the  current  study,  they  are  of  some  potential  importance. 
Eigenvalue  solutions  {see  Reference  32,  for  example)  should  contain 
a  Sturm  sequence  check  in  this  case.  This  matter  is  further 
discussed  in  Section  7. 


3.8  DESIGN  SCALING 

The  process  for  determining  the  scale  factor  A  in  equation 
(3.68)  is  trivial  in  the  case  of  stress  and  deflection  constraints. 
Stresses  and  deflections  are  simply  inversely  proportional  to  the 
member  sizes.  For  constraints  on  natural  frequencies,  the  situation 
is  not  as  simple.  This  can  be  seen  by  observing  that  when  the  mass 
and  stiffness  matrices  are  scaled  in  the  Rayleigh  quotient 


}T [K]  {xi> 
{x;.}TIM]  {xi> 


(3.69) 


the  frequency  will  remain  unchanged.  This  problem  was  noted  by 
Kuisalaas  in  Reference  33. 

Frequencies  can  be  modified,  however,  by  scaling  a  portion 
of  the  mass  and  stiffness  matrices.  Assume  that  the  global  mass 
and  stiffness  matrices  have  the  forms 
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(3.70a) 


[M] 

[K] 


m 

[M°]  + 

m 

>V  + 


(3.70b) 


where  the  [Ml  and  [K  ]  matrices  correspond  to  mass  and  stiffnesses 
which  remain  fixed  during  the  design  scaling  and  optimization  process. 
The  a^'s  are  a  set  of  relative  design  variables  which  are  to  be 
scaled  by  A  to  achieve  a  required  natural  frequency,  w  . 

Introduction  of  the  unknown  scale  factor  creates  a  new 
eigenvalue  problem  with  A  as  the  eigenvalue.  Rearrangement  of 
equation  (3.69)  using  equations  (3.70)  and  w.  =  w  ,  {x.}  =  {xe} 

1  K  1  K 

gives 


A [A] {xR}  =  [B] {xR} 

(3.71a) 

m  2 

where 

(A]  =  l  a  (w  [M  ]  -  [K  ]  ) 
i=l  1  K  1  1 

(3.71b) 

2 

and 

[BJ  =  [K  J  "  w_[M  ] 

O  R  O 

(3.71c) 

The  smallest  real  eigenvalue,  A*,  of  the  above  problem  (if  one 

exists)  will  be  the  design  scale  factor  which  places  the  design 

exactly  on  the  constraint  boundary,  w.  -  w_  =  0 . 

1  K 

To  demonstrate  the  above  scaling  procedure,  the  optimal 
design  of  a  tapered  bar  with  a  tip  mass  was  studied.  Figure  3,4 
shows  both  a  continuous  model  and  a  discrete  model  of  the  bar. 
Reference  34  gives  the  optimal  mass  distribution  as 


M  (x) 


M_  cosh' 


(SjL)/ 


cosh" 


(31x) 


(3.72) 


where  Mq  is  the  tip  mass,  L  is  the  length. 


with 


P 
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Figure  3.4.  Tapered  Bar  with  Tip  Mass. 


being  the  mass  density  and  E  is  Young's  modulus.  The  fundamental 
frequency,  is  fixed.  The  total  optimal  mass  of  the  continuous 

bar  is 


MT  =  Mo  (3.73) 

which  is  the  integral  of  equation  (3.72),  For  illustrative 
pumoses,  the  following  data:  L  =  90  in.,  Mq  =  .05  lb. -sec. 2/in. , 

P  =  2.59  x  10"4  lb. -sec. 2/in. 4 ,  E  =  10.3  x  106  psi,  w_  -  2578  rad/sec 

”  o 

are  used.  Use  of  these  values  gives  MT  =  .104  lb. -sec.  /in.  The 
discrete  model  of  Figure  3.4(b)  consists  of  two  axial  bar  elements 
of  equal  length. 

Using  the  above  material  properties,  constraint  curves  can 
be  generated  for  various  choices  of  the  fundamental  frequency,  u>R. 
Figure  3.5  shows  frequency  contours  in  the  M1  -  M2  design  space 
for  various  choices  of  the  eigenvalue  **  pw2£.2y  6E.  The  =  .056 
curve  corresponds  to  the  Reference  34  value  of  w  .  For  convex 

R 

feasible  regions  such  as  those  indicated  by  Figure  3.5,  the  scaling 
of  the  design  along  a  vector  from  the  origin  to  a  feasible  point 
is  a  well  defined  process. 

Figure  3.6  shows  contours  of  constant  eigenvalues  in  the 
second  mode.  In  this  case  the  optimal  design  will  be  controlled 
by  side  constraints  on  M1  and  M2  as  well  as  the  eigenvalue 
requirement.  Thus,  trends  that  are  noted  for  the  fundamental  mode 
may  not  carry  over  for  higher  modes. 

The  value  of  the  fixed  mass,  Mq,  also  strongly  affects  the 

frequency  constraint.  In  Figure  3.7  three  values  of  M  were 

o 

selected  (.05,  ,10,  and  .20)  and  the  contours  of  A^  =  .056  were 
determined.  For  =  .05  the  optimal  design  is  -  .076, 

M2  =  .0  34  and  hence  M^,  =  .110  which  is  about  six  percent  higher 
than  the  continuous  case.  As  Mq  approaches  zero,  the  optimal 
design  problem  becomes  meaningless  without  side  constraints. 
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Figure  3.5.  Mode  1  Eigenvalue  Contours  for  Two  Element  Bar 
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Figure  3.6,  Mode  2  Eigenvalue  Contours  for  Two  Element  Bar. 


056  (fixed) 
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Figure  3.7.  Effect  of  Tip  Mass  Variations  on  Mode  1  Eigenvalue 
for  Two  Element  Bar . 


The  following  design  strategy  achieves  minimum  weight  for 
a  given  fundamental  frequency  requirement; 

1.  Form  the  characteristic  determinant  for  the  system 
shown  in  equation  (3.71),  treating  A  as  the  eigenvalue. 

2.  Find  the  smallest  positive  root  A*  and  the  corresponding 
eigenvector  {X  *}. 

K 

3.  Compute  the  weight  using  A*  as  the  scale  factor. 

4.  Compute  the  Lagrangian  densities,  e^,  for  each  element. 

5.  Resize  according  to  M.  =  M.  ,  ,  /\  e .  t 

^  i  new  i  old  /  1  i( 

6.  Normalize  the  Mi  ngw  design  and  go  to  Step  2. 

7.  Repeat  Steps  2  -  5  a  specified  number  of  times  or  until 
the  weight  increases. 


The  above  scaling  and  resizing  procedure  tends  to  equalize 
the  energy  {Lagrangian  density)  in  all  the  elements  and  at  the 
same  time  brings  the  design  to  the  required  fundamental  frequency. 

It  can  be  shown  that  the  Lagrangian  densities  are  directly  proportional 
to  the  rates  of  change  of  the  eigenvalues  with  respect  to  the  element 
names.  Thus  uniform  Lagrangian  density  implies  equal  sensitivity 
to  design  changes.  The  energies  can  sometimes  be  negative  for  a 
given  element.  This  implies  that  increasing  the  element  causes 
a  reduction  in  the  frequency.  Note  that  the  element  resizing 
formula  uses  the  absolute  value  of  the  element  energies. 

In  the  bar  example,  the  mass  and  stiffness  matrices  are 
tridiagonal.  The  Sturm  sequence  of  polynomials 


P1  =  1 
P2  =  C1 


Pi  ci-lPi-l  "  bi-2  Pi-2 


i=2, . . . ,n+l 


are  generated  where  the  c^ ' s  are  the  diagonal  and  the  b^ 1 s  are 
the  upper  and  lower  diagonals.  Sign  alternations  need  to  be 
checked  to  assure  that  the  smallest  eigenvalue  is  obtained.  The 
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root  is  determined  by  using  the  golden  section  minimization 
35 

technique  on  the  absolute  value  of  the  polynomial  p 
Reference  32  contains  a  listing  of  the  eigenroutine.  This  soft¬ 
ware  package  also  contains  a  Givens-Householder36  routine  to 
reduce  equations  (3.71)  to  tridiagonal  form  for  more  complex 
problems  which  are  not  tridiagonal.  The  golden  section  technique 
together  with  Sturm  sequence  checks  is  more  efficient  than  the 
usual  bisection  or  secant  procedure. 

Figure  3.8  shows  scaling  and  resizing  steps  in  design  space 
for  the  two  element  bar .  One  scaling  took  the  design  from  point  A 
to  point  B  which  lies  on  the  constraint.  An  energy  resizing 
scaling  step  then  brought  the  design  to  point  C.  Figure  3.9  shows 
how  resizing  proceeds  from  a  different  starting  point.  Additional 
resizings  and  scalings  from  point  D  tend  to  create  a  "zig-zag" 
path  in  which  little  additional  weight  is  removed. 

The  above  design  scaling  technique  for  frequency  requirements 
needs  further  testing  before  it  is  included  in  the  ADDRESS  program. 
The  results  in  this  section  are  intended  to  show  that  frequency 
requirements  can  be  met  in  a  noniterative  manner.  The  combination 
of  scaling  and  energy  resizing  appears  to  have  considerable 
promise  for  frequency  constraints. 
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Figure  3.8.  Frequency  Constraint  Resizing  of  Two  Element  Bar 
(One  Energy  Resizing) . 
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Figure  3.9.  Frequency  Constraint  Resizing  of  Two  Element  Bar 
(Two  Energy  Resizing) . 


SECTION  4 

ALGORITHM  FOR  DAMAGE  REQUIREMENTS 

In  this  section  the  algorithm  for  the  optimal  design  of 
damage  tolerant  structures  is  presented.  After  the  overview  of 
the  approach,  the  organization  of  the  calculations  is  presented. 

This  is  followed  by  a  discussion  of  the  use  of  reanalysis  procedures 
and  a  description  of  the  computer  code  which  implements  the 
procedure . 

4.1  OVERVIEW  OF  THE  APPROACH 

A  conventional  approach  to  optimization  and  damage  assess¬ 
ment  may  proceed  along  the  following  steps: 

(1)  Define  the  undamaged  structure. 

(2)  Solve  the  optimized  problem  - 

find  {d} 

such  that  W({d})  +  minimum 

and  G.{d}>^0  i  =  1,...,N 

1  —  o 

where  {d}  =  design  variable  vector 

W  =  merit  function  (weight) 
t  h 

G^  =  i  constraint  on  undamaged  structure 

N  =  number  of  constraints, 
c 

(3)  Perform  vulnerability  analysis  of  damage  configuration  using 
the  {d}  determined  in  (2). 

(4)  Adjust  {d}  to  obtain  satisfactory  performance  of  the  damaged 
structure . 

The  last  step  could  involve  engineering  judgment  and  the  use  of 
criteria  based  on  experience  with  similar  design  problems. 

The  more  general  approach  taken  herein  is  to  include 
additional  constraints  for  the  damage  configuration  in  the  problem 
statement  as  follows: 

(1)  Define  the  undamaged  and  damaged  configurations. 

(2)  Solve  the  optimization  problem  - 
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find  d 

such  that  W(d)  -*■  minimum 

and  G. (d)  >0  i  *  1,...,N 

1  —  ""  c 

and  g-^td)  _>  0  i  =  1,...,N  ^ 

m 

g  .  (d)  >  0  i  =  1 . N 

mi  -  —  cm 

th  tti 

where  -  i  constraint  on  the  k  damage  model 
m  =  number  of  damage  models 
N£,jc  =  number  of  constraints  on  the  kth  damage  model 

and  d,  W,  G. ,  N  are  as  before.  The  constraint  conditions  include 
-  1  c 

requirements  on  stress,  deflection,  natural  frequency,  and  limits 
on  member  sizes.  As  shown  in  Figure  4.1,  this  general  approach 

also  includes  the  possibility  for  modification  to  the  applied 
loads  in  the  damage  models. 

4.2  ORGANIZATION  OF  THE  CALCULATIONS 

Calculations  to  perform  the  optimization  for  damage  tolerance 
requirements  are  arranged  as  given  in  the  flow  description  below 
for  a  typical  design  step: 

(1)  Assemble  stiffness  (and  mass,  if  required)  for  the  latest 
design. 

(2)  Perform  static  solutions  for  all  design  conditions  not 
containing  damage.  Save  the  factored  stiffness  matrix  for 
use  in  (4) . 

(3)  If  required,  perform  solutions  for  natural  frequencies  and 
mode  shapes  of  the  undamaged  structure. 

(4)  For  each  design  condition  containing  damage,  define  the 
change  in  stiffness  and  solve  for  the  response  of  the  damaged 
structures  using  reanalysis  techniques.  Find  changes  in 
frequencies  and  mode  shapes  of  the  damaged  structures,  if 
required. 
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Figure  4.1.  Damage  Optimization  Concept. 


(5) 


Compute  element  forces,  stress,  etc.  of  the  undamaged  and 
damaged  configurations  in  preparation  for  resizing.  If  an 
element  is  damaged  in  a  certain  design  condition,  generate 
null  stress  vector  (or  other  appropriate  responses)  for  that 
condition • 

(6)  Perform  resizing,  element  by  element,  using  optimality 

criteria . 

(7)  Go  to  (1)  to  begin  the  next  redesign  cycle,  if  required. 

The  above  computation  scheme  integrates  damage  reanalysis 
into  the  optimization  procedure  in  such  a  way  that  the  presence  of 
damage  in  certain  design  conditions  is  transparent  to  the  resizing 
segment  of  the  program.  Solutions  for  damage  conditions  are 
returned  as  additional  response  vectors,  for  which  noncritical 
stress  states  are  to  be  generated  in  any  damaged  elements.  Thus, 
the  optimization  can  proceed  as  usual  with  behavioral  constraints 
for  damaged  design  conditions  treated  in  an  analogous  manner  as 
for  additional  loading  cases  in  an  undamaged  structure. 

4.3  THE  ADDRESS  PROGRAM  STRUCTURE 

The  computer  program  ADDRESS,  developed  under  this  project, 
was  designed  from  the  outset  to  be  easily  modified  and  to 
incorporate  state  of  the  art  numerical  methods.  The  code  was 
developed  in  standard  ANSI  FORTRAN  and  contains  no  features  which 
preclude  its  portability.  The  code  was  written  in  a  structured 
form  and  arranged  in  modular  routines  with  liberal  use  made  of 
internal  documentation. 

From  a  conceptual  point  of  view,  the  overall  program  control 
organization  is  shown  in  Figure  4.2.  The  "Executive  Control",  a 
main  program,  calls  the  "Analytic/Optimization  Control"  which,  in 
turn,  performs  the  analysis  and  optimization  by  calls  to  the 
energy  (stress) ,  displacement,  and  frequency  control  segments. 

The  Executive  Control  is  detailed  in  Figure  4.3.  Its  primary 
functions  are  to  set  up  the  structural  model  (both  original  and 
damaged)  and  to  control  the  optimization  process. 
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Figure  4.2.  Overall  Control  Organization  of  the  ADDRESS  Program. 


EXECUTIVE  CONTROL 


—  INPUT/INITIALIZATION 

—  ELEMENT  STIFFNESS,  MASS  AND 

CHANGE  IN  STIFFNESS,  MASS 

—  MATRIX  ASSEMBLY  AND 

BOUNDARY  CONDITIONS 

1 —  DESIGN  STRATEGY  {STEP  BY  STEP 
DECISIONS) 

—  ANALYSIS/OPTIMIZATION  CONTROL 


Figure  4.3.  ADDRESS  Executive  Control. 
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The  design  strategy,  shown  in  Figure  4.4,  currently  consists 
of  a  preselected  number  of  cycles  of  energy  mode  resizing,  followed 
by  displacement  mode  resizing  if  required.  Frequencies  are  checked 
at  each  step.  Stress  and  deflection  scaling  takes  place  to  assume 
a  feasible  design  after  each  step. 

Details  of  the  energy,  displacement,  and  frequency  mode 
controls  are  given  in  Figures  4.5,  4.6,  and  4.7,  respectively.  The 
energy  and  displacement  modes  are  nearly  identical  to  the  strategy 
of  the  OPTSTATCOMP  program.  The  frequency  mode  control  involves 
the  use  of  the  SIMIT  technique  for  the  undamaged  structure  and  the 
vibration  reanalysis  procedure  of  Section  3.5  herein.  Further 
details  of  these  procedures  are  given  in  the  ADDRESS  User's 
Manual . 
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ANALYSIS/ 

OPTIMIZATION 

CONTROL 


STRESS  (ENERGY  MODE) 

DISPLACEMENT  CONSTRAINT 
MODE 

FREQUENCY  MODE 

DESIGN  VECTOR  AVERAGING 
(STEP  SIZE  MODIFICATION) 


Figure  4.4,  ADDRESS  Analysis/Opt imization  Control. 
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ENERGY  MODE 


ANALYSIS 


EQUATION  SOLUTION 
(UNDAMAGED  CASES) 

ITERATIVE  REANALYSIS 
(DAMAGED  CASES) 

ELEMENT  FORCE,  STRESS 
ENERGY  (ALL  CASES) 


t 


SCALING  FOR  STRESS 
CONSTRAINT  (CRITICAL  CASE) 

ELEMENT  RESIZING 


Figure  4.5.  ADDRESS  Energy  Mode  Control. 
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displacement 

MODE 


ANALYSIS 

h 


“  EQUATION  SOLUTION 
(UNDAMAGED  CASES) 

—  ITERATIVE  REANALYSIS 
(DAMAGED  CASES) 


DUMMY  LOAD  CASES 
( DAMAGED/UNDAMAGED ) 


ELEMENT  RESIZING 


Figure  4.S.  ADDRESS  Displacement  Mode  Control. 
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FREQUENCY 

MODE 


ANALYSIS 


DECOMPOSITION 
{UNDAMAGED  CASE) 

SIMULT.  ITERATION 
(UNDAMAGED  CASE) 

FREQUENCY/MODE 
RE ANALYSIS 
'DAMAGED  CASES) 


Figure  4.7.  ADDRESS  Frequency  Mode  Control. 
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SECTION  5 


APPLICATIONS  OF  THE  ADDRESS  PROGRAM  TO 
REPRESENTATIVE  STRUCTURES 

The  ADDRESS  program  was  applied  to  several  representative 
structures.  Analysis,  reanalysis,  and  optimization  results  are 
presented  in  this  section.  Comparisons  are  made  with  programs 
obtained  from  the  Air  Force  such  as  DANALYZ  (vibration  program) 
and  OPTSTATCOMP  (optimization  program  with  composites  capability) . 
The  effects  of  damage  on  the  structural  response  and  the  optimal 
member  sizes  are  noted.  Results  are  presented  for  the 
following  structures: 

•  Ten  bar  truss 

•  Simple  wing  box 

•  Intermediate  complexity  wing  (ICW) 

A  major  design  and  analysis  study  of  the  A- 7  outer  wing  panel 
was  also  conducted.  These  results  are  given  in  the  next  section. 

5 . 1  TRUSS  RESULTS 

The  optimal  design  of  the  ten  bar  truss  shown  in  Figure  5.1 
has  been  studied  in  detail  by  a  number  of  investigators  (see 
Reference  30,  for  example).  Two  independent  load  cases  (LC1,  LC2 ) 
shown  in  Figure  5.2  are  applied  together  with  the  two  damage 
cases  (DC1 ,  DC2 )  shown  in  Figure  5.3. 

The  constraints  are  that  member  stresses  not  exceed  20,000  psi 
and  that  nodal  deflections  be  less  than  or  equal  to  2.0  inches. 

In  the  first  damage  case,  95%  of  the  material  is  removed  from 
member  9,  and  in  the  second  case  95%  of  the  stiffness  in  members 
2,  5,  and  6  is  removed.  The  material  properties  are  based  on 
aluminum. 

Six  truss  designs  were  studied  as  shown  in  Table  5.1.  The 
first  truss  design  has  all  the  member  sizes  equal  and  is  the 
lowest  weight  uniform  design  which  satisfies  all  the  constraints 
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Figure  5.1.  Ten-Bar  Truss  Geometry. 
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Figure  5.2.  Ten-Bar  Truss  Load  Cases. 
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Figure  5.3.  Ten-Bar  Truss  Damage  Cases. 


TABLE  5.1 

DESCRIPTION  OF  TRUSS  DESIGNS 


Truss 

Desiqn 

Description 

Weight 

(lb) 

1 

Uniform 

8267 

2 

Optimal  for  LC1 

5077 

3 

Optimal  for  LC2 

4823 

4 

Optimal  for  LC1,  2 

5455 

5 

Optimal  for  LC1,  DC1 

6866 

6 

Optimal  for  LC1,  DC1,  2 

9581 

TABLE  5.2 

SUMMARY  OF  TRUSS  RESULTS 

<5  =  Maximum  Deflection 

max 


NOE  -  Number  of  Overstressed  Elements 


6 

6 

max 

max 

TD 

LC 

DC 

NOE 

(in) 

TD 

LC 

DC 

NOE 

(in) 

1 

1 

1 

0 

2.5 

4 

1 

1 

5 

5.6 

2 

0 

2.3 

2 

0 

2.1 

2 

1 

0 

2.6 

2 

1 

6 

6.2 

2 

1 

2.2 

2 

4 

-9.4 

2 

1 

1 

5 

10.8 

5 

1 

1 

0 

2.0 

2 

1 

2.0 

2 

7 

15.2 

2 

1 

6 

11.6 

2 

1 

3 

-10.4 

2 

5 

-283.5 

2 

7 

14.3 

3 

1 

1 

3 

8.9 

6 

1  , 

1 

0 

1.9 

2 

0 

2.4 

2 

0 

1.7 

2 

1 

5 

8.9 

2 

1 

0 

2.0 

2 

4 

-7.8 

2 

1 

1.6 

65 


under  the  first  load  case.  Designs  2,  3,  and  4  axe  optimal  for 
the  load  cases  indicated  without  damage  to  any  of  the  members. 

The  fifth  and  sixth  designs  include  the  indicated  damage  conditions 
when  the  structure  is  loaded  under  the  first  condition.  As  can  be 
seen  from  the  table,  the  addition  of  damage  considerations  causes 
an  increase  in  weight,  just  as  multiple  load  cases  increase  the 
weight  of  the  optimal  design. 

Table  5.2  is  a  summary  of  load  and  damage  tolerance  studies 
performed  on  the  six  truss  designs.  The  uniform  design  performs 
quite  well  compared  to  the  optimal  designs  which  do  not  include 
damage  tolerance.  The  designs  which  are  optimal  for  only  one  load 
case  are  quite  sensitive  to  changes  in  loading  or  addition  of 
damage.  The  sixth  design,  which  is  optimal  for  the  first  load 
case  and  both  damage  conditions,  is  also  satisfactory  under  the 
second  load  case  when  only  the  first  damage  is  imposed. 

The  above  calculations  were  performed  with  the  aid  of  a 
small  scale  computer  program  which  is  documented  in  Appendix  E  of 
the  ADDRESS  User's  Manual.  Its  primary  purpose  is  to  reanalyze 
the  damaged  truss  structure  using  simple  iteration.  Insight  was 
gained  into  the  number  of  iterations  required  to  obtain  convergence. 
In  general,  no  more  than  five  iterations  are  required  to  achieve 
stress  and  displacements  within  5%  of  those  obtained  from  an  exact 
analysis.  This  applies  to  damages  of  the  level  used  in  this 
problem.  For  damages  which  introduce  a  kinematic  instability  in 
the  truss,  the  iterations  diverge. 

5.2  SIMPLE  WING  BOX 

To  check  the  efficiency  of  the  ADDRESS  calculations,  a 
small-scale  finite  element  model  of  a  straight,  rectangular  wing 
box  was  developed.  The  model  is  shown  in  Figure  5.4  and  is 
composed  of  24  bar,  shear  panel,  and  membrane  cover  elements. 

There  are  a  total  of  48  unrestrained  degrees  of  freedom  in  the 
model  corresponding  to  the  three  degrees  of  freedom  at  each  of  the 
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Figure  5.4.  Simple  Wing  Box  Finite  Element  Model. 


sixteen  unrestrained  nodes.  The  covers  in  bay  3  were  reduced  in 
mass  by  25%  and  in  stiffness  by  50%  to  simulate  a  ballistic  damage. 

The  wing  was  subjected  to  loads  of  10,000  lb.  at  the  four 
tip  nodes  in  the  negative  z-direction.  Constraints  consisted  of 
the  requirement  that  no  member  have  stresses  in  excess  of 
20,000  psi  and  that  displacements  be  no  greater  than  2.0  inches. 

The  structure  is  designed  for  minimum  weight  such  that  the  stress 
and  deflection  requirements  are  not  exceeded  for  the  original  or 
the  damaged  wing  box. 

Three  cycles  of  energy  resizing  were  employed  with  scaling 
to  meet  the  constraints.  The  element  weights  are  shown  in 
Table  5.3.  Elements  were  linked  such  that  top  and  bottom  covers 
are  of  the  same  thickness,  shear  panels  vary  together,  and  bars 
(posts)  are  all  the  same.  The  results  show  that  after  one 
iteration  the  weight  has  reached  a  near  constant  value  of  307  lb. 

Table  5.4  gives  comparisons  between  using  a  complete 
analysis  at  each  step  and  using  the  reanalysis  techniques 
discussed  in  Section  3.  An  accuracy  requirement  of  0.005  on  both 
the  displacements  and  the  mode  shapes  was  required.  The  table 
shows  that  although  the  computing  times  are  small,  the  computing 
time  is  reduced  by  75%  for  the  static  analysis  and  by  over  50%  for 
the  vibration  analysis. 

5.3  INTERMEDIATE  COMPLEXITY  WING 

The  intermediate  complexity  wing  model  of  Reference  30  is 
shown  in  Figure  5.5.  It  contains  a  total  of  158  membrane,  bar, 
and  shear  panel  elements.  To  test  the  convergence  of  the  static 
and  vibration  reanalysis  iterations,  25  elements  were  damaged  in 
the  leading  edge  tip  area  as  shown  in  the  figure.  Two  separate 
load  conditions  corresponding  to  a  subsonic,  forward  center  of 
pressure  load  and  a  supersonic,  near  uniform  pressure  loading 
were  assumed.  The  maximum  deflection  occurred  at  the  indicated 
deflection  control  point  for  the  supersonic  loading. 
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TABLE  5 . 3 


SIMPLE  WING  BOX  WEIGHTS  FOR  ENERGY  AND 
DISPLACEMENT  RESIZING 


Element 

Number 

Type 

Connectivity 

Initial 

Weight 

Cycle  1 
Weight 

Cycle  2 
Weight 

Cycle  3 
Weight 

MA 

MB 

MC 

MD 

1 

C 

i 

3 

7 

5 

37. 

.088 

8. 

,177 

7. 

.537 

7 

.505 

2 

C 

2 

4 

8 

6 

8. 

,177 

7. 

.537 

7 

.505 

3 

C 

5 

7 

11 

9 

30. 

,542 

30. 

,935 

30 

.979 

4 

c 

6 

8 

12 

10 

30. 

,542 

30. 

.935 

30. 

.  979 

5 

C 

9 

11 

15 

13 

37. 

,719 

37. 

.493 

37 

.463 

6 

C 

10 

12 

16 

14 

37. 

,719 

37. 

.493 

37 

.463 

7 

c 

13 

15 

19 

17 

50. 

,250 

50. 

.634 

50 

.650 

8 

c 

14 

16 

20 

18 

50. 

,250 

50. 

.634 

50 

.653 

9 

SP 

1 

2 

6 

5 

14. 

.835 

6. 

,730 

6. 

,  745 

6 

.745 

10 

SP 

3 

4 

8 

7 

11 

SP 

5 

6 

10 

9 

12 

SP 

7 

8 

12 

11 

13 

SP 

9 

10 

14 

13 

14 

SP 

11 

12 

16 

15 

15 

SP 

13 

14 

18 

17 

16 

SP 

15 

16 

20 

19 

17 

B 

1 

2 

0 

0 

.  593 

014 

* 

,010 

.  010 

18 

B 

3 

4 

0 

0 

19 

B 

5 

6 

0 

0 

20 

B 

7 

8 

0 

0 

21 

B 

9 

10 

0 

0 

22 

B 

11 

12 

0 

0 

23 

B 

13 

14 

0 

0 

24 

B 

15 

16 

0 

0 

1 

1 

' 

Total  Weights 

(all  elements) 

420. 

128 

307. 

328 

307. 

.238 

307, 

.234 

C  =  Cover 
SP  =  Shear  Panel 
B  =  Bar 
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TABLE  5.4 


COMPARISONS  OF  ANALYSIS  AND  REANALYSIS 
COMPUTER  TIMES  FOR  SIMPLE  WING  BOX 


Computing  Time  ( sec ) 

Calculation 

Initial 

Cycle  1 

Cycle  2 

Cycle  3 

Total 

Static  Analysis  (sec) 

.  31 

.  30 

.30 

.29 

1.20 

Static  Reanalysis  (sec) 

.08 

.07 

.07 

.06 

.28 

Eigen  Analysis  (sec) 

.21 

.24 

.24 

.25 

.94 

Eigen  Reanalysis  (sec) 

.11 

.10 

.10 

.09 

.40 
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Table  5.5  shows  the  iteration  history  and  compares  it  to 
the  "exact"  analysis  of  the  damaged  structure.  Figure  5.6  shows 
the  error  in  the  tip  deflection  as  a  function  of  computing  time. 

The  computing  cost  for  five  iterations  is  seen  to  be  about 
one-half  that  of  the  "exact"  analysis.  At  the  same  time,  the  tip 
deflection  determined  by  the  iteration  is  within  1%  of  the  exact 
value.  The  time  savings  on  larger  problems  is  even  greater  since 
iteration  involves  only  successive  forward  and  back  substitutions. 

To  check  the  vibration  analysis  capability,  the  ICW  model 
was  run  using  the  Air  Force's  DAN AL Y Z  program  and  compared  to  the 
simultaneous  iteration  (SIMIT)  technique  used  in  ADDRESS.  Several 
damages  were  imposed  on  the  wing  covers  as  shown  in  Figure  5.7. 
Table  5.6  summarizes  the  frequency  calculations  using  DANALYZ  and 
SIMIT  to  compute  the  first  four  modal  frequencies.  Good 
correlation  was  obtained  in  all  cases.  Mode  shapes  compared  to 
three  decimal  places.  Mode  shapes  appear  to  be  relatively 
insensitive  to  structural  damage  for  this  problem  (see  Figures  5.8 
through  5.19) , 

Table  5.7  gives  results  of  the  iteration  procedure  for  the 
ICW  damage  condition  of  Figure  5.5.  The  iteration  column  of 
Table  5.7  indicates  the  number  of  times  which  modes  and  frequencies 
have  been  recomputed  using  the  iteration  scheme  of  Section  3. 

It  is  interesting  to  note  that  not  all  frequencies  converge  at  the 
same  rate.  For  mode  one  and  three,  one  iteration  was  sufficient  to 
achieve  a  high  degree  of  accuracy.  Mode  two,  however,  required 
three  iterations  to  achieve  similar  accuracy.  The  efficiency  of 
the  reanalysis  is  also  quite  good.  The  table  shows  that  computing 
time  can  be  reduced  by  over  50%  for  this  problem. 

5.4  COMPARISONS  WITH  OPTSTATCOHP 

The  efficiency  and  accuracy  of  the  ADDRESS  program  was 
compared  with  that  of  OPTSTATCOMP  for  several  of  the  problems  in 
the  previous  subsections.  In  some  cases,  direct  comparisons  are 
not  possible  since  different  minimum  sizes  or  number  of  resizings 
were  used. 
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TABLE  5.5 


CONVERGENCE  OF  STATIC  REANALYSIS  FOR  ICW 
TIP  DEFLECTION  AT  DEFLECTION  CONTROL  POINT 


Iteration 

Number 

Tip  Deflection 
(inches) 

undamaged 

20.3512 

1 

21.1567 

5 

21.6716 

10 

21.6911 

15 

21.6916 

20 

21.6916 

exact* 

21.7018 

*" Exact"  in  the  sense  that  a  complete  analysis 
was  performed  without  iteration. 
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Figure  5,6.  Static  Reanalysis  Computing  Times  for  ICW. 
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Figure  5.7.  Damage  Zones  on  the  Intermediate  Complexity  Wing  (ICW) . 
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Figure  5.8.  First  Bending  Mode  for  ICW  (No  Damage). 
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Figure  5.9.  Second  Bending  Mode  for  ICW  (No  Damage) . 
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Figure  5.10.  First  Torsion  Mode  for  ICW  (No  Damage) 


80 


Figure  5. IX.  Third  Bending  Mode  for  ICW  (No  Damage). 
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Figure  5.12.  First  Bending  Mode  for  ICW  (Zone  1  Damage). 
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Figure  5.13.  Second  Bending  Mode  for  ICW  (Zone  1  Damage). 
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Figure  5.14.  First  Torsion  Mode  for  ICW  (Zone  1  Damage) 
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Figure  5.15.  Third  Bending  Mode  for  ICW  (Zone  1  Damage). 
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Figure  5.16.  First  Bending  Mode  for  ICW  (Four  Damages). 
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Figure  5.17.  Second  Bending  Mode  for  ICW  (Four  Damages). 
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Figure  5.18.  First  Torsion  Mode  for  ICW  (Four  Damages). 
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Figure  5.19.  Third  Bending  Mode  for  ICW  {Four  Damages). 


TABLE  5.7 


FREQUENCY  REANALYSIS  RESULTS  FOR  ICW 


Iteration 

CP  Sec. 
(CYBER  74) 

fl 

(Hz) 

f2 

(Hz) 

f  3 
(Hz) 

0* 

- 

40.93 

132.76 

168.63 

1 

2.42 

45.94 

132.58 

189.03 

2 

3.58 

45.70 

131.15 

189.80 

3 

4.74 

45.67 

130.97 

189.59 

4 

5.90 

45.66 

130.91 

189.56 

exact** 

12.17 

45.66 

130.87 

189.52 

*Frequencies  at  iteration  zero  are  for  the 
undamaged  wing. 

**Exact  frequencies  and  CP  times  for  a  complete 
eigen  analysis  of  the  damaged  structure. 
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Table  5.8  shows  the  results  for  the  ten-bar  truss  example. 

The  ADDRESS  program  and  the  OPTSTATCOMP  program  were  both  run 
for  five  cycles  of  energy  resizing.  Reference  30  shows  that  if 
more  than  ten  iterations  are  used,  the  weight  should  converge  to 
about  5062  pounds  when  minimum  sizes  of  0.1  square  inches  are  used 
for  this  load  case. 

Tables  5.9  and  5.10  show  the  ADDRESS  and  OPTSTATCOMP  results, 
respectively,  for  the  metal  ICW  model  discussed  earlier.  The 
stress  optimization  using  both  programs  yields  a  final  design 
weight  after  five  cycles  of  approximately  193  pounds.  The 
ADDRESS  program  seems  to  produce  a  slightly  lower  weight  than 
OPTSTATCOMP  after  one  cycle.  Direct  stress  allowables  are  taken 
to  be  25  ksi  in  both  tension  and  compression. 

A  stress  and  deflection  optimization  comparison  is  shown  in 
Tables  5.11  and  5.12  for  the  ICW.  The  first  four  cycles  are 
based  on  energy,  and  the  remaining  ones  are  resizing  for  the 
deflection  requirement.  In  the  ADDRESS  program,  if  the  weight 
more  than  doubles  in  the  displacement  resizing  mode,  then  execution 
is  terminated  and  the  previous  design  is  taken  to  be  the  final  one. 
The  OPTSTATCOMP  program  continued  to  resize  for  five  cycles  of 
displacement  calculations,  even  though  the  doubling  of  weight 
criterion  was  met.  It  then  returned  to  the  energy  based  design  of 
cycle  5. 

A  composite  ICW  Model  based  on  the  data  of  Reference  19  was 
also  run.  Only  stress  requirements  were  included.  Tables  5.13  and 
5.14  show  the  comparative  results.  Due  to  the  higher  minimum 
layer  requirement  used  in  ADDRESS  (.0053)  compared  to  OPTSTATCOMP 
(.005),  the  weights  of  the  cover  elements  shown  in  Table  5.13  are 
higher.  As  with  the  metal  ICW,  the  composite  version  was  resized 
for  two  load  conditions. 

Figure  5.20  shows  the  total  number  of  layers  in  the  top  and 
bottom  covers  obtained  from  the  two  programs.  In  Figure  2.21  the 
number  of  layers  is  broken  down  into  numbers  in  the  0°,  90°,  ±45° 
directions.  The  0°  direction  is  taken  to  be  the  direction  of  the 
center  spar. 
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TABLE  5 . 8 


COMPARISON  OF  ADDRESS  AND  OPTSTATCOMP 
TEN-BAR  TRUSS  DESIGNS 


Weight  in  Pounds 


Element 

Number 

ADDRESS 

Weight 

OPTSTATCOMP 

Weight 

1 

1106.4 

1093.2  ‘ 

2 

3.6 

17.6 

3 

861.7 

884.1 

4 

530.5 

524.8 

5 

3.6 

.4 

6 

3.6 

24.8 

7 

434.4 

466.5 

8 

1067.0 

1049.8 

9 

1061.1 

1049.7 

10 

5.1 

35.2 

Total  Weight 

5077.0 

5146.1 
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TABLE  5.9 


ADDRESS  RESULTS  FOR  METAL  ICW  WITH  STRESS  CONSTRAINT* 


Weight  in  Pounds 


Element 

Type 

Initial 

Weight 

Cycle  1 
Weight 

Cycle  2 
Weight 

Cycle  3 
Weight 

Cycle  4 
Weight 

Cycle  5 
Weight 

Covers 

328.98 

191.33 

189.24 

182.29 

177.66 

171.35 

Shear  Panels 

102.44 

23.62 

23.09 

22.54 

22.29 

21.78 

Bars 

6.74 

.25 

.24 

.22 

.22 

.21 

Total  Weight 
(all  elements) 
- 1 

438.16 

215.20 

212.57 

205.05 

200.17 

193.34 

TABLE  5.10 

OPTSTATCOMP  RESULTS  FOR  METAL  ICW  WITH  STRESS  CONSTRAINT* 


Weight  in  Pounds 


Element 

Type 

Initial 

Weight 

Cycle  1 
Weight 

Cycle  2 
Weight 

Cycle  3 
Weight 

Cycle  4 
Weight 

Cycle  5 
Weight 

Covers 

328.98 

193.98 

189.04 

182.20 

177.61 

171.23 

Shear  Panels 

102.44 

23.94 

23.04 

22.51 

22.26 

21.75 

Bars 

6.74 

.26 

.24 

.  22 

.22 

.21 

Total  Weight 
(all  elements) 

438.16 

218.18 

212.32 

204.93 

200.09 

193.19 

*Allowables  set  at  25  ksi. 
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TABLE  5.11 

ADDRESS  RESULTS  FOR  METAL  ICW  WITH  STRESS*  AND  DEFLECTION**  CONSTRAINTS 
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*Allowables  set  at  60  ksi. 

**Deflections  limited  to  4.0  inches  in  z-direction. 


TABLE  5.13 


ADDRESS  RESULTS  FOR  COMPOSITE  ICW 
WITH  STRESS  CONSTRAINT 


Element 

Initial 

Cycle  1 

Cycle  2 

Cycle  3 

Cycle  4 

Type 

Weight 

Weight 

Weight 

Weight 

Weight 

Covers 

124.94 

117.13 

98.15 

9  3.27 

92.89 

Shear  Panels 

70.73 

43.42 

27.62 

27.06 

27.38 

Bars 

4.65 

2.75 

1.49 

1.45 

1.47 

Total  Weight 
(all  elements) 

200.33 

163.30 

127.26 

121.78 

121.74 

Weight  in  Pounds 


TABLE  5.14 

OPTSTATCOMP  RESULTS  FOR  COMPOSITE  ICW 
WITH  STRESS  CONSTRAINT 


Element 

Initial 

■Cycle  1 

Cycle  2 

Cycle  3 

Cycle  4 

Type 

Weight 

Weight 

Weight 

Weight 

Weight 

Covers 

124.94 

86.03 

76.03 

71.24 

71.24 

Shear  Panels 

70.73 

33.80 

28.42 

27.10 

27.59 

Bars 

4.65 

2.02 

1.53 

1.43 

1.47 

Total  Weight 
(all  elements) 

200.33 

121.85 

105.98 

99.77 

100.30 

Weight  in  Pounds 


94 


NOTE: 


One  OPTSTATCOMF  layer 
=  .005  in. 


One  ADDRESS  layer 
=  .0053  in. 


Top  figures  are  from 
ADDRESS:  bottom 
figures  {  )  are  from 
OPTSTATCOMP. 


Figure  5.20.  Total  Number  of  Layers  in  the  Top  Skin  of 
Composite  ICW. 
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NOTE: 


Top  figures  are  from 
ADDRESS;  bottom 
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Figure  5,21 


Distribution  of  Fibers  in  0°,  90°,  ±45° 
Directions  in  the  Top  Skin  of  Composite  ICW . 


Table  5.15  summarizes  the  execution  times  for  the  two  programs 
on  the  problems  just  discussed.  One  major  reason  for  the  faster 
ADDRESS  time  is  in  the  coding  of  the  decomposition  of  the  stiffness 
matrix.  Subroutine  GAUSS  was  replaced  with  subroutine  LLT  which 
has  the  DO- loops  reorganized.  This  section  of  code  in  ADDRESS 
now  runs  about  four  times  as  fast  as  the  OPTSTATCOMP  code.  If 
reanalysis  had  been  used  in  these  examples,  the  running  times 
would  have  been  improved  even  more. 
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TABLE  5.15 

EXECUTION  TIMES  FOR  OPTSTATCOMP  AND  ADDRESS 


Problem 

Execution  Times  -  CYBER  175 

OPSTATCOMP 

ADDRESS 

TEN  BAR  TRUSS 

Stress  and  Displacement 
Optimization 

0.9 

0.5 

METAL  WING 

Stress  Optimization 

38.9 

23.8 

METAL  WING 

Stress  and  Displacement 
Optimization 

53.5 

28.7 

COMPOSITE  WING 

Stress  Optimization 

76.9 

51.7 
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SECTION  6 

DESIGN  STUDIES  FOR  A-7D  OUTER  WING  PANEL 

To  demonstrate  the  ADDRESS  computer  program,  a  major  task 
of  this  work  was  to  perform  a  design  vulnerability  and  optimization 
study  on  an  existing  aircraft  lifting  surface  component.  This 
section  describes  the  UDRI  and  Vought  Corporation  efforts  to 
model,  resize,  and  further  refine  the  design  of  an  A-7D  outer  wing 
panel.  Both  metal  and  composite  models  are  considered. 

The  work  described  in  this  section  involved  several  distinct 
subtasks.  Vought,  acting  as  a  subcontractor  to  the  UDRI, 
participated  as  a  consultant  in  the  selection  of  the  A7D  outer  wing 
panel  as  the  candidate  structure  for  redesign.  During  the  analysis 
task  Vought  supplied  UDRI  with  geometry,  loads,  material  properties, 
and  ballistic  damage  requirements.  The  UDRI  then  set  up  the  finite 
element  model  and  performed  optimization  studies  using  ADDRESS. 
Member  sizes  were  then  supplied  to  Vought  and  additional  design 
criteria  were  applied.  The  final  designs  were  checked  again 
using  the  ADDRESS  program  to  obtain  final  weights  and  to  check  the 
stress  levels.  A  test  plan  and  design  drawings  were  generated  by 
Vought  for  fabrication  and  testing  purposes.  These  items  are 
documented  as  a  separate  contract  item. 

6.1  FINITE  ELEMENT  MODELING  OF  THE  A-7D  OUTER  WING  PANEL 

The  A-7  aircraft  is  an  attack  fighter  designed  for  carrier- 
based  operations.  Figure  6.1  shows  the  overall  aircraft  arrange¬ 
ment.  The  portion  of  the  wing  structure  which  was  modeled  in  this 
study  is  shaded  in  the  figure. 

Figure  6,2  shows  the  A-7  wing  structure  in  some  detail.  The 
portion  of  the  wing  inboard  of  the  hinge  line  is  all  metal.  The 
outer  wing  panel  can  be  replaced  and  is  of  metal  or  composite 
construction.  Further  structural  detail  of  the  outer  panel  is 
shown  in  Figure  6.3.  It  is  this  general  arrangement  of  Figure  6.3 
which  formed  the  basis  of  the  finite  element  model. 
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The  wing  structure  was  tested  originally  for  the  static 
loads  described  in  Table  6.1.  A  review  of  these  design  loads 
determined  that  the  condition  LAD-1  produces  the  highest  skin 
stresses  and  was  therefore  included  as  a  load  case.  The  second 
loading  used  was  the  RP0-2A  condition  which  is  primarily  a  torsion 
loading. 

Figures  6.4  and  6.5  show  the  shear,  moment,  and  torque 
diagrams  for  the  LAD— 1  and  the  RP0-2A  load  cases,  respectively. 

These  distributions  apply  to  the  entire  wing  from  the  fuselage  to 
the  tip.  Only  the  portions  of  the  curves  outboard  of  the  wing 
fold  line  (yw  =  142  in.)  were  used  to  generate  the  discrete  loads. 

An  overview  of  the  UDRI  finite  element  model  is  shown  in 
Figure  6.6.  The  seven-spar  arrangement  was  used  in  both  the  metal 
and  composite  design  studies.  The  model  contains  113  nodes  and 
453  finite  elements.  Nodes  104  through  113  are  fully  restrained 
to  simulate  a  fixed  condition  at  the  wing  hinge  line.  The  nodes 
numbered  in  the  figure  apply  to  the  top  surface  with  the  missing 
numbers  applying  to  the  lower  surface.  The  wing  tip  tapers  to  a 
line  of  nodes  1  through  7,  and  hence  there  are  no  corresponding 
lower  surfaces  at  span  station  93 . 

Two  damage  zones  are  also  indicated  in  Figure  6.6.  Damage 
zone  number  one  consists  of  the  eight  covers  (upper  and  lower) 
shaded  near  the  tip.  The  second  damage  zone  is  near  the  root  and 
consists  of  four  covers.  Two  other  cases  were  also  considered  for 
several  computer  runs.  In  one  of  these,  rib  and  spar  elements  were 
also  damaged  as  well  as  the  indicated  covers.  The  final  damage 
case  consisted  of  reducing  the  stiffness  of  all  the  covers  from 
root  to  tip  between  the  third  and  fifth  spar  from  the  leading  edge. 
In  all  damage  cases,  75%  of  the  element  thickness  was  removed  to 
simulate  the  effect  of  ballistic  damage  to  the  zone. 

Table  6.2  shows  the  nodal  coordinates  and  the  three  components 
of  the  applied  loads  for  both  load  cases.  The  point  loads  were 
iteratively  determined  in  such  a  way  as  to  match  as  closely  as 
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Figure  6.4.  A- 7  Load  Case  Number 
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29.73" 
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113  Nodes  (104  through  113  are 

fully 

restrained) 

• 

309  Degrees  of  freedom 

• 

1-94  Membrane  covers 

• 

95-142  Shear  panels  (ribs) 

• 

143-194  Shear  panels  (spars) 

# 

195-296  Bars  (rib  caps) 

• 

297-400  Bars  (spar  caps) 

• 

401-453  Bars  (posts) 
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Figure  6.6.  Overview  of  A- 7  Outer  Wing  Panel 
Finite  Element  Model. 
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TABLE  6 . 2 

NODAL  AND  APPLIED  LOAD  DATA  FOR  A-7D  OUTER  WING  PANEL 


Node 

No. 

X 

in 

T7 

in 

Z 

in 

F  - 
X 

LCl 

F  - 

X 

LC2 

F  - 

y 

LCl 

v: 

LC2 

F  -  LCl 
z 

F  -  LC2 
z 

I 

74.  60 

92.70 

0.00 

0.00 

0*00 

0* 

00 

0- 

00 

17.83 

15.26 

19.77 

16.92 

2 

80.89 

93.23 

0*00 

3 

84.04 

93.31 

0.00 

21*13 

18*08 

4 

87.19 

93.39 

0.00 

22.29 

19.08 

5 

89.61 

93*37 

0.00 

23*07 

19.74 

6 

92*02 

93*  34 

0.00 

20.54 

17,58 

7 

104.53 

92.99 

0,00 

16,4  7 

14 , 10 

B 

72.00 

89.51 

0*43 

160 .47 

158*48 

9 

72,23 

09*34 

-1*76  | 

160,47 

158*48 

10 

78.07 

89.35 

1,39 

178*12 

175,89 

11 

76.07 

89, 35 

-1,39 

178.12 

175*89  1 

12 

81-54 

89.35 

1.46 

190*47 

182.85 

13 

01.54 

89.35 

-1.46 

190,47 

182*85 

14 

85.00 

89. 35 

1.54 

202.83 

200*30 

15 

85.00 

B9 . 35 

-1*54 

202*83 

200,30 

16 

87*71 

89. 34 

1.55 

195*77 

193*31 

17 

87*71 

89.34 

-1.55 

195*77 

193,32 

IB 

90.42 

89*33 

1.56 

218*71 

215.99 

19 

90.  42 

89.33 

-1,56 

218*71 

215-99 

20 

102.97 

89.30 

1*13 

259.30 

256.06 

21 

102.97 

89*30 

-1*13 

259*30 

256,06 

22 

63.26 

78*64 

*60 

528,12 

474,33 

23 

63.52 

78.44 

-1*95 

528*12 

474,33 

24 

CD 

75,18 

1*61 

475.12 

421.24 

25 

67.78 

75.18 

-1*62 

475.12 

421,24 

26 

71.73 

72,67 

1.75 

427*43 

373*45 

27 

71*73 

73,67 

-1*75 

427*43 

373,45 

23 

75.67  . 

72*16 

1*88 

474 . 44 

320*35 

29 

75,67 

72*16 

-1.88 

474*44 

320*35 

30 

79,09 

70,00 

1.88 

342,64 

298*49 

31 

79.09 

70,80 

-1,88 

342.64 

288,49 

32 

82.50 

69.44 

1*87 

300,25 

246.02 

33 

82,50 

69*44 

-1*87 

300*25 

246*02 

34 

92,84 

65*30 

1*31 

178.37 

123,89 

35 

92,84 

65*30 

-1.31 

178,37 

123.99 

36 

49*71 

61.79 

.74 

265*72 

775*63 

37 

49.99 

61,56 

-2*33 

265*72 

775.63 

3fl 

53-91 

56*10 

1*94 

362*88 

667*49 

39 

53.91 

56,10 

-1*96 

362,88 

667.49 

40 

59.27 

1  53*96 

2,12 

468,13 

571,16 

41 

59*27 

53*  96 

-2*13 

468*13 

571.16 

42 

64.63 

51,  B1 

2*29 

565*29 

473*01 

43 

64.63 

51.81 

-2*29 

565.29 

473,01 

44 

69.31 

4  9,94 

2*23 

|  638.15 

399*40 

45 

69.31 

49*94 

-2*24 

638.15 

399.40 

46 

73*99 

48.07 

2, IB 

727,21 

309,44 

47 

73*99 

48,07 

-2*18 

V 

727.21 

309*44 

48 

83.92 

44*10 

1*58 

889.14 

109*86 
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TABLE  6.2  (Continued) 


Node 

Ho. 

X 

In 

Y 

in 

2 

In 

Fx 

-  LC1 

Fx  -  LC2 

F 

y 

-  LCl 

F 

y 

-  3X2 

F  -  LCl 
z 

F  -  LC2 

Z 

49 

8  3.92 

4  4.10 

-1*50 

0. 

00 

0, 

00 

0 

.00 

0, 

00 

839,14 

109.86 

50 

29*80 

37*04 

*  95 

-201,38 

1213.24 

51 

30.19 

36.75 

-2.88 

-201.38 

1213.24 

52 

36*27 

31*83 

2,37 

149.15 

1039,66 

S3 

36.27 

31.83 

-2.  38 

149.15 

1038,66 

54 

42.63 

27.41 

2*61 

499,67 

864.09 

55 

42.63 

27.41 

-2.62 

499,67 

864*09 

56 

48.98 

22*98 

2.06 

850,20 

699*51 

57 

48.98 

22,98 

-2.86 

850.20 

689,51 

56 

55*47 

20,43 

2.72 

1135,00 

547.67 

59 

55.47 

20.43 

-2.72 

1135*00 

547.67 

60 

61.96 

17.87 

2.57 

1463*62 

384,00 

61 

61.96 

17.07 

-2.57 

1463.62 

384.00 

62 

71.25 

14.16 

1.97 

1857.96 

187,61 

63 

71*25 

14.16 

-1.97 

1857.96 

187*61 

64 

o*oc 

0.00 

1.22 

-126.60 

ir8,95 

65 

*  56 

“•  37 

-3.  84 

-126.60 

268,95 

66 

17.85 

6.48 

2.82 

-25.14 

214.19 

67 

17*85  , 

6.48 

-2.83 

-25.14 

214*19 

60 

28*98 

6.48 

3.00 

63.65 

166*23 

69 

2  8.98 

6.48 

-3,05 

63,65 

166,23 

70 

40.11 

6,48 

3.19 

143,97 

122.93 

71 

40,11 

6.48 

-3.26 

143.97 

122.93 

72 

48.  76 

6.48 

2.96 

207.38 

88*70 

73 

48,76 

6,48 

-2*99 

207,38 

88,70 

74 

57.40 

6,48 

2.72 

270.80 

54,4  7 

75 

57.40 

6.48 

-2,  72 

270.80 

54.47 

76 

66.29 

2.39 

2.33 

346.89 

13.40 

77 

66.29 

2,39 

-2.33 

346.89 

13.40 

78 

13,12 

-  -  02 

2,67 

-29,36 

216,48 

79 

13.12 

- .  02 

-2.91 

-29.36 

216.40 

80 

16,54 

2.60 

2,84 

-16.68 

209-98 

81 

16.54 

,2,60  , 

-2.80 

-16.68 

209,98 

82 

27.86 

2,46 

3*07 

59*42 

169,56 

83 

27.86 

2.46 

-3*04 

59.42 

168.56 

84 

39.18 

2,  32 

3.29 

152.43 

118*36 

85 

39,18 

2*32 

-3.27 

152.43 

118*36 

86 

47,01 

2.32 

3.00 

207. 38 

89.70 

87 

47.01 

2.  32 

-3*04 

207. 38 

99*70 

88 

54.83 

2.32 

2.71 

266*57 

56.76 

89 

54.83 

2.32 

-2,80 

266*57 

56.76 

90 

16.54 

,58 

2*87 

-8.22 

205.07 

91 

16.54 

,58 

-2.84 

* 

-8.22 

205.07 

92 

27,86 

.58 

3.08 

* 

76.33 

159.43 

93 

27, B6 

,58 

-3,05 

76,33 

159.43 

94 

39.18 

.SB 

3*29 

160-88 

113,00 

95 

39*18 

.58 

-3.27 

160,89 

113,80 
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TABLE  6 . 2  (Continued) 


Node 
No « 

X 

in 

Y 

in 

t 

in 

F 

X 

LCl 

FX  - 

*  LC2 

F  - 

y 

■  LCl 

F  ■ 

Y 

-  LC2 

F  - 
z 

LCl 

F  - 
z 

LC2  i 

96 

46,23 

.58 

3.03 

0. 

00 

0.00 

0.00 

0.00 

211 

.61 

86 

.42 

97 

46.23 

.58 

-2.98 

211 

.61 

86 

.42 

98 

53.27 

.58 

2,77 

262 

.34 

59 

.04 

99 

53.27 

.58 

-2.6B 

262 

.  34 

59 

.04 

100 

59.  37 

.47 

2.62 

304 

.62 

36 

.22 

101 

59.37 

,47 

-2.62 

304 

.62 

36 

.22 

102 

65.45 

.39 

2.36 

34  6 

.89 

13 

.40 

103 

65.45 

.  39 

-2.36 

346 

.89 

13 

.40 

104 

16.54 

-1,12 

2.70 

0 

.00 

0 

.00 

105 

16.54 

-1.12 

-2.68 

106 

27.86 

-1.12 

2.99 

107 

27.86 

-1.12 

-2.98 

108 

39. 18 

-1.12 

3.28 

109 

39.18 

-1.12 

-3.27 

110 

45.45 

-1.12 

3.00 

111 

45.45 

-1.12 

-2.98 

112 

51,71 

-1.12 

2.72 

1 

f 

113 

51.71 

-1.12 

-2.70 

1 

F  ,  F  ,  F  loads  are  in  pounds, 
x  y  z 

LC1  =  Load  Case  1 

LC2  =  Load  Case  2 
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possible  the  shear,  moment,  and  torque  diagrams.  The  cumulative 
effect  of  the  loads  at  the  root  matches  the  load  diagram  exactly; 
however,  due  to  the  skewed  pattern  of  the  nodes,  the  load  data 
in  the  midsection  of  the  wing  matches  only  to  within  +10%. 

Figures  6.7  through  6.10  provide  further  details  of  the 
model.  The  wing  taper  and  tip  modeling  discussed  above  can  be 
clearly  seen  in  Figure  6.7.  Figure  6.8  shows  the  cover  numbering 
scheme.  The  rib  and  spar  panel  numbers  are  shown  in  Figures  6.9 
and  6.10,  respectively. 

6.2  MATERIAL  PROPERTIES  AND  DAMAGE  CASES 

Material  properties  for  the  metal  A7D  outer  panel  were  based 
on  typical  7075-T6  aluminum: 

E  =  Young's  modulus  =  10.3  x  10®  psi 
v  =  Poisson's  ratio  =  .33 

3 

p  =  Mass  density  =  .1  lbf/in 
G  =  Shear  modulus  =  3.9  x  10®  psi 

with  the  following  stress  allowables: 

aTU  =  Allowable  stress  in  tension  =  72  ksi 
aCY  =  Allowable  stress  in  compression  =  63  ksi 
aCR  =  Allowable  stress  in  shear  =  43  ksi. 

The  composite  design  was  based  on  5208/T-300  graphite/epoxy  tape 
with  a  layup  consisting  of  50%  (0°  plies),  40%  (±45°  plies),  and 
10%  (90°  plies) : 

=  12.2  x  10®  psi 
E 2  =  4.5  x  10®  psi 
V  =  .4 

p  =  .057  lbf/in3 
G  =  2.7  x  10®  psi 

with  the  following  stress  allowables: 

°TUl  =  102,0  ksi 
cr ipy 2  —  15.0  ksi 
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Figure  6.7.  A- 7  Outer  Wing  Panel  Mode  Numbers.  (Sheet 
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Figure  6.7.  A- 7  Outer  Wing  Panel  Mode  Numbers.  (Sheet  2  of  2) 
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Figure  6.8.  A- 7  Outer  Wing  Panel  Cover  Elements.  (Sheet  1  of  2) 
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Figure  6.8.  A- 7  Outer  Wing  Panel  Cover  Elements.  (Sheet  2  of  2) 
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Figure  6.9.  A- 7  Outer  Wing  Panel  Rib  Panel  Elements. 
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Figure  6.10.  A- 7  Outer  Wing  Panel  Spar  Panel  Elements. 


°CY1  =  45 ' 6  kSi 
°CY2  *  45'6  k3i 

tCr  ■  30’5  ksi 

where  the  1  and  2  subscripts  denote  properties  in  the  longitudinal 
and  transverse  directions,  respectively. 

It  was  found  that  significant  factors  in  the  optimal  design 
were  the  selection  of  the  allowable  skin  stress  level  and  the 
maneuver  load  factors  after  damage.  The  allowable  stresses  in 
tension  for  the  skin  after  damage  were  computed  in  terms  of  the 
allowable,  cTU,  of  the  undamaged  structure  from  the  expression 


TU 


=  «d(W>0TD 


(6.1) 


where  6,  is  a  damage  reduction  factor  based  on  impact  fracture  or 
a  37 

residual  strength  considerations.  In  developing  the  specific 
criteria.  References  38  and  39  were  used  extensively. 

For  the  metal  design,  the  o*y  for  the  cover  elements  was 
computed  by  equation  (6.1)  using  <5^  =  .16,  Nz  =  7 ,  and  N*  =  3. 

Two  sets  of  allowables  were  used  for  the  composite  designs.  The 
nominal  set  was 


a  * 

TUI 

°TU2 

°CYl 

aCY2 

T  * 


21.4  ksi 
3.2  ksi 
45.6  ksi 
45.6  ksi 

30.5  ksi 


A  second  set 

°TU1 

aTU2 

aCYl 

aCY2 

"T  *  — 


or  "reduced" 

11.3  ksi 

4.3  ksi 

45.6  ksi 

45.6  ksi 
27.0  ksi 


set  consisted  of 


Details  of  the  selection  of  these  values  and  their  relationship  to 
ballistic  damage  are  given  in  Reference  37. 
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Table  6.3  lists  the  damage  cases  studied.  Cases  1,  2,  3,  and 
8  are  for  the  metal  wing.  In  Case  1  the  wing  was  sized  for  stress 
constraints  considering  two  applied  loads,  and  no  ballistic  damage. 

Case  2  included  damage  to  the  two  cover  areas  indicated  in 
Figure  6.6.  In  Case  3  the  ribs  and  spars  were  also  damaged.  The 
last  metal  design.  Case  8,  allowed  for  damage  to  all  cover  elements 

*  from  the  root  to  the  tip  between  the  third  and  fifth  spar. 

Cases  4  through  7  represent  studies  for  the  composite  panel 
design.  Case  4  included  no  damage,  Case  5  included  the  cover 
damage  of  Case  2,  Case  6  was  the  same  as  Case  5  but  used  the  above 
''reduced"  allowables,  and  Case  7  included  the  same  damage  as  Case  3. 

6.3  DESIGN  RESULTS 

The  ADDRESS  computer  program  was  used  to  resize  the  metal  and 
composite  A7D  panels.  Table  6.4  summarizes  the  weights  for  the 
eight  cases  discussed  in  the  previous  section.  Three  cycles  of 
energy  resizing  were  used  and  the  weights  computed  after  each  cycle. 

It  is  noted  in  Table  4  that  the  initial  weights  are  not  the  same 
in  all  cases.  This  is  due  to  the  fact  that  the  designs  were 
uniformly  scaled  to  reach  the  stress  allowables  for  the  various 
cases . 

Figures  6.11,  6.12,  and  6.13  give  an  indication  for  the 
Case  2  metal  design  how  the  material  was  redistributed  among  the 
covers,  spars,  and  ribs.  The  various  levels  of  shading  correspond 

# 

to  the  indicated  range  of  the  member  sizes. 

.  Figure  6 . 14  summarizes  the  Case  5  composite  design  of  the 

covers.  Ribs  and  spar  element  sizes  were  very  similar  to  the  metal 
design. 

Using  criteria  for  local  failure  set  forth  in  Reference  37, 
the  Case  2  and  Case  5  optimal  designs  were  further  modified.  The 
values  of  the  member  sizes  so  obtained  are  given  in  Table  6.5. 

The  weights  shown  in  Table  6.6  summarize  the  "idealized”  (ADDRESS) 
and  "design"  (Reference  37)  weights.  The  first  set  of  idealized 
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TABLE  6.3 

A7D  OPTIMIZATION  CASE  DESCRIPTIONS 


Case 

Description 

1 

Metal,  no  damage 

2 

Metal,  2  damages  (covers  only) 

3 

Metal,  2  damages  (all  elements) 

4 

Composite,  no  damage 

5 

Composite,  2  damages  (covers  only) 

6 

Case  5  with  reduced  allowables 

7 

Case  6,  2  damages  (all  elements) 

8 

Metal,  1  damage  in  second  and 
third  bay  from  root  to  tip 
(covers  only) 
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TABLE  6.4 


A-7D  CASE  WEIGHTS  FOR  ENERGY  RESIZINGS 


Case 

Elements 

Initial 

Cycle  1 

Cycle  2 

Cycle  3 

Weight 

Weight 

Weight 

Weight 

1 

C 

325.1 

88.3 

73.4 

69.7 

SP 

181.5 

36.7 

28.0 

25.9 

B 

92.9 

19.3 

15.4 

14.6 

Total 

599.5 

144.3 

116.8 

110.2 

2 

c 

350.1 

114.9 

94.1 

82.7 

SP 

195.5 

46.5 

33.6 

28.9 

B 

100.1 

25.5 

20.5 

19.1 

Total 

645.7 

186.9 

148.2 

130.7 

3 

c 

361.1 

125.2 

107.3 

99.2 

SP 

201.6 

50.8 

38.0 

33.9 

B 

103.2 

27.2 

21.6 

20.0 

Total 

665.9 

203.2 

166.9 

153.1 

4 

C 

254.6 

60.6 

53.6 

52.7 

SP 

248.4 

31.8 

24.9 

23.7 

B 

127.2 

17.2 

13.4 

11.8 

Total 

630.2 

109.6 

91.9 

88.2 

5 

C 

277.4 

84.3 

75.1 

70.9 

SP 

270.7 

43.2 

31.2 

27.8 

B 

138.6 

24.3 

18.6 

16.0 

Total 

686.7 

151.8 

124.9 

114.7 

6 

C 

529.2 

152.3 

138.7 

131.7 

SP 

516.4 

74.2 

52.3 

46.9 

B 

264.4 

40.0 

29.9 

25.3 

Total 

1310.0 

266.5 

220.9 

203.9 

7 

C 

546.7 

165.6 

155.6 

151.8 

SP 

533.5 

81.4 

59.2 

54.3 

B 

273.2 

42.6 

30.9 

25.2 

Total 

1353.4 

289.6 

245.7 

231.3 

8 

C 

350.4 

111.4 

92.5 

82.4 

SP 

195.6 

44.5 

32.7 

28.4 

B 

100.2 

24.6 

20.4 

19.1 

Total 

646.2 

180.5 

145.6 

129.9 

C  =  Covers 
SP  =  Shear  Panels 
B  =  Bars 
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Figure  6.11.  A-7D  Covers  Resized  for  Damage  (Metal). 
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Figure  6.12.  A-7D  Spars  Resized  for  Damage  (Metal). 


* 
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Figure  6.13.  A-7D  Rib  Resized  for  Damage  (Metal). 
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Figure  6.14.  A-7D  Covers  Resized  for  Damage  (Composite). 


TABLE  6.5 


FINAL  DESIGNS  FOR  METAL  AND  COMPOSITE  A-7D 
OUTER  WING  PANEL 


Element 

No. 

Type 

Metal 

Design 

Composite 

Design 

Element 

No. 

Type 

Metal 

Design 

Composite 

Design 

1 

4 

.432 

.530 

47 

4 

.069 

.110 

2 

.343 

.435 

48 

.080 

.125 

3 

.052 

.065 

49 

.040 

.060 

4 

.051 

.  060 

50 

.054 

.085 

5 

.  052 

.  065 

51 

.247 

.263 

6 

.051 

.060 

52 

.  304 

.259 

7 

.  052 

.  065 

53 

.253 

.280 

8 

.051 

.060 

54 

.301 

.268 

9 

.052 

.  065 

55 

.269 

.299 

10 

.051 

.  060 

56 

.  319 

.279 

11 

.121 

.170 

57 

.272 

.315 

12 

.121 

.170 

58 

.346 

.284 

13 

.036 

.040 

59 

.  124 

.  180 

14 

.036 

.040 

60 

.141 

.205 

15 

.059 

.068 

61 

.039 

.  055 

16 

.055 

.060 

62 

.036 

.050 

17 

.059 

.068 

63 

.096 

.170 

18 

.055 

.060 

64 

.096 

.170 

19 

.  059 

.068 

65 

.138 

.220 

20 

.  054 

.060 

66 

.  127 

.200 

21 

.059 

.068 

67 

.284 

.  420 

22 

.054 

.060 

68 

.273 

.  400 

23 

.  036 

.040 

69 

.  340 

.510 

24 

.036 

.040 

70 

.  365 

.545 

25 

.  045 

.070 

71 

.309 

.  455 

26 

.  036 

.060 

72 

.  327 

.  475 

27 

.  076 

.103 

73 

.167 

.205 

28 

.  059 

.080 

74 

.073 

.130 

29 

.  080 

.  098 

75 

.090 

.170 

30 

.  059 

.  080 

76 

.099 

.180 

31 

.078 

.  098 

77 

.157 

.250 

32 

.  059 

.  080 

78 

.145 

.225 

33 

.078 

.098 

79 

.274 

.  390 

34 

.059 

.  080 

80 

.286 

.420 

35 

.036 

.040 

81 

.554 

.  825 

36 

.  036 

.  050 

82 

.531 

.785 

37 

.072 

.105 

83 

.504 

.745 

38 

.  082 

.125 

84 

.563 

.830 

39 

.134 

.049 

85 

.036 

.050 

40 

.  097 

.143 

86 

.036 

.050 

41 

.134 

.157 

87 

.134 

.225 

42 

.097 

.145 

88 

.117 

.215 

43 

.  132 

.157 

89 

.173 

.270 

44 

.096 

.145 

90 

.163 

.  255 

45 

.130 

.159 

91 

.227 

.335 

46 

1 

1 

.  094 

.147 

92 

► 

.220 

.325 
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TABLE  6.5  (Continued) 


Element 

No. 

Type 

Metal 

Design 

Composite 

Design 

93 

4 

.  886 

1.335 

94 

* 

.912 

1.375 

95 

5 

.  098 

.082 

96 

.  040 

.  040 

97 

.040 

.  040 

98 

.040 

.  040 

99 

.040 

.040 

100 

.036 

.035 

101 

.  042 

.039 

102 

.050 

.050 

103 

.050 

.050 

104 

.050 

.  050 

105 

.050 

.  050 

106 

.036 

.035 

107 

.036 

.035 

108 

.063 

.063 

109 

.063 

.063 

110 

.063 

.063 

111 

.063 

.063 

112 

.036 

.035 

113 

.074 

.  074 

114 

.090 

.090 

115 

.090 

.090 

116 

.090 

.090 

117 

.090 

.090 

118 

.112 

.112 

119 

.036 

.  035 

120 

.100 

.100 

121 

.100 

.100 

122 

.100 

.100 

123 

.  100 

.100 

124 

.036 

.035 

125 

.036 

.035 

126 

.036 

.035 

127 

.036 

.035 

12  8 

.036 

.035 

129 

.036 

.035 

130 

.036 

.  035 

131 

.036 

.035 

132 

.036 

.035 

133 

.036 

.035 

134 

.036 

.035 

135 

.036 

.035 

136 

.036 

.035 

137 

.036 

.  035 

138 

.055 

.048 

139 

.  039 

.035 

140 

i 

f 

.036 

.035 

Element 

No. 

Type 

Metal 

Design 

Composite 

Desiqn 

141 

5 

.036 

.035 

142 

.036 

.035 

143 

.  801 

.722 

144 

.063 

.063 

145 

.063 

.063 

146 

.063 

.063 

147 

.063 

.063 

148 

.063 

.063 

149 

.043 

.041 

150 

.049 

.047 

151 

.063 

.  063 

152 

.063 

.063 

153 

.063 

.063 

154 

.063 

.063 

155 

.063 

.063 

156 

.036 

.035 

157 

.036 

.035 

158 

.063 

.063 

159 

.063 

.063 

160 

.063 

.063 

161 

.063 

.063 

162 

.063 

.063 

163 

.  036 

.035 

164 

.072 

.061 

165 

.080 

.  080 

166 

.  080 

.080 

167 

.080 

.080 

168 

.080 

.080 

169 

.080 

.  080 

170 

.075 

.063 

171 

.037 

.035 

172 

.125 

.125 

173 

.125 

.125 

174 

.125 

.125 

175 

.125 

.125 

176 

.125 

.125 

177 

.060 

.054 

178 

.041 

.035 

179 

.055 

.035 

180 

.077 

.062 

181 

.137 

.054 

182 

.055 

.077 

183 

.092 

.123 

184 

.134 

.127 

185 

.169 

.162 

186 

.174 

.167 

187 

.227 

.210 

188 

1 

f 

.312 

.292 
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TABLE  6.5  (Continued) 


Element 

No. 

Type 

Metal 

Design 

Composite 

Design 

Element 

No. 

n 

Type 

1  Metal 
Design 

Composite 

Design 

189 

5 

.338 

.314 

305 

2 

.  036 

.035 

190 

.223 

.228 

306 

.036 

.035 

191 

.338 

.  346 

307 

.060 

.  060 

192 

.347 

.352 

308 

.  060 

.059 

193 

.036 

.035 

309 

.036 

.035 

194 

.063 

.  055 

310 

.036 

.035 

195 

2 

.036 

.  035 

311 

.  036 

.035 

196 

.  036 

.0  35 

312 

.036 

.035 

197 

.  036 

.0  35 

313 

.036 

.0  35 

198 

.036 

.035 

314 

.036 

.035 

199 

.  036 

.035 

315 

.036 

.035 

200 

.036 

.035 

316 

.036 

.035 

201 

.036 

.035 

317 

.036 

.035 

202 

.  046 

.040 

318 

.036 

.035 

203 

.036 

.035 

319 

.036 

.035 

1 

320 

.036 

.035 

321 

.036 

.035 

322 

.036 

.035 

274 

.036 

.035 

323 

.036 

.035 

275 

.  060 

.06  7 

324 

.036 

.035 

276 

.036 

.035 

325 

.036 

.  035 

277 

.036 

.035 

326 

.043 

.035 

278 

.  026 

.035 

327 

.041 

.035 

279 

.  036 

.  035 

328 

.036 

.035 

280 

.036 

.  035 

329 

.039 

.035 

281 

.036 

.035 

330 

.  036 

.035 

282 

.  036 

.  035 

331 

.037 

.035 

283 

.036 

.035 

332 

.036 

.035 

284 

.036 

.035 

333 

.036 

.035 

285 

.036 

.035 

334 

.036 

.035 

286 

.036 

.  035 

335 

.036 

.035 

287 

.036 

.035 

336 

.036 

.035 

288 

.  036 

.035 

337 

.036 

.035 

289 

.036 

.  035 

338 

.036 

.0  35 

290 

.036 

.  035 

339 

.036 

.  035 

291 

.036 

.035 

340 

.107 

.073 

292 

.036 

.0  35 

341 

.096 

.068 

293 

.036 

.035 

342 

.063 

.044 

294 

.036 

.035 

343 

.  099 

.068 

295 

.  036 

.  035 

344 

.  083 

.  057 

296 

.036 

.035 

345 

.115 

.080 

297 

.734 

.659 

346 

.109 

.076 

298 

.  422 

.386 

347 

.102 

.073 

299 

.036 

.035 

348 

.103 

.074 

300 

.  036 

.035 

349 

.  096 

.  072 

301 

.036 

.0  35 

350 

.104 

.077 

302 

.036 

.  035 

351 

.045 

.  035 

303 

.036 

.035 

352 

.059 

.  040 

304 

1 

' 

.  036 

.  035 

353 

1 

> 

.036 

.035 
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TABLE  6.5  (Continued) 


Element 

No. 

Type 

Metal 

Design 

Composite 

Design 

Element 

No. 

Type 

Metal 

Design 

' 

Composite 

Design 

35  4 

2 

.  036 

.  035 

405 

2 

.036 

.035 

355 

.066 

.047 

406 

.036 

.035 

356 

.045 

.  035 

407 

.036 

.  035 

357 

.189 

.145 

408 

.036 

.035 

358 

.167 

.129 

409 

.036 

.035 

359 

.686 

.555 

410 

,036 

.035 

360 

.654 

.529 

411 

.036 

.035 

361 

.460 

.338 

412 

.036 

.035 

362 

.470 

.344 

413 

.036 

.035 

363 

.195 

.133 

414 

.036 

.035 

364 

.212 

.144 

415 

.036 

.035 

365 

.036 

.035 

416 

.  036 

.035 

366 

.036 

.  035 

417 

.036 

.035 

367 

.054 

.050 

418 

.036 

.035 

368 

.073 

.  066 

419 

.036 

.035 

369 

.042 

.047 

420 

.036 

.035 

370 

.057 

.  059 

421 

.036 

.035 

371 

.050 

.  053 

422 

,036 

.035 

372 

.094 

.098 

423 

.036 

.035 

373 

.124 

.106 

424 

.036 

.035 

374 

.118 

.097 

425 

.036 

.035 

375 

.121 

.106 

426 

.036 

.0  35 

376 

.084 

.  075 

427 

.036 

.035 

377 

.068 

.063 

428 

.036 

.  035 

378 

.052 

.048 

429 

.036 

.0  35 

379 

.176 

.127 

430 

.045 

.040 

380 

.159 

.  114 

431 

.036 

.035 

381 

.14  8 

.118 

4  32 

.036 

.035 

382 

.154 

.126 

433 

.036 

.035 

383 

.237 

.190 

434 

.036 

.035 

384 

.217 

.172 

435 

.036 

.035 

385 

.397 

.284 

4  37 

.036 

.035 

386 

.358 

.251 

438 

.036 

.035 

387 

.281 

.194 

4  39 

.036 

.035 

388 

.378 

.277 

440 

.127 

.130 

389 

.149 

.099 

441 

.205 

.212 

390 

.125 

.077 

442 

.036 

.035 

391 

.298 

.199 

443 

.036 

.035 

392 

.  313 

.207 

444 

.036 

.035 

393 

.  473 

.337 

445 

.142 

.145 

394 

.527 

.375 

446 

.248 

.253 

395 

1.206 

.969 

447 

.036 

.035 

396 

1.279 

1.016 

448 

.036 

.035 

397 

.  036 

.035 

449 

.036 

.  035 

398 

.036 

.035 

450 

.036 

.035 

399 

.036 

.035 

451 

.036 

.  035 

400 

.036 

.  035 

452 

1 

,036 

.035 

401 

.036 

.035 

453 

.036 

.035 

402 

.036 

.035 

403 

.036 

.035 

404 

.036 

.035 
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TABLE  6.6 

A7D  WEIGHT  SUMMARY 


Case 

Metal 

Wingbox 

■ 

Composite 

Wingbox 

Idealized  Weight 
{static  criteria  but  with 
no  ballistic  damage) 

116.9  lbs 

88.2  lbs 

Idealized  Weight 
(UDRI  results) 

130.7  lbs 

114.7  lbs 

Pinal  Design  Weight 
(Vought  results) 

160.6  lbs 

123.3  lbs 

* 


130 


weights  represents  the  undamaged  structure  which  is  sized  according 
to  static  stress  requirements.  The  second  set  of  idealized  weights 
represents  the  structure  sized  according  to  ballistic  damage 
criteria  {only  covers  damaged) .  The  idealized  weights  are  based 
solely  upon  the  ADDRESS  program  output  for  the  undamaged  and  damaged 
wingboxes.  The  design  weights  reflect  a  mix  of  UDRI  output  and 
Vought  baseline  wingbox  data  as  reflected  on  the  design  drawings. 

The  design  weights  for  the  metal  and  the  composite  wingboxes 
increased  due  to  the  method  by  which  the  structure  was  sized.  The 
method  used  made  an  element-by-element  comparison  of  Vought  input 
data  with  ADDRESS  output.  The  larger  of  the  two  values  was  then 
used  in  the  design.  This  method  assured  that  the  wingboxes  would 
meet  all  static,  fatigue  and  fracture  requirements  as  well  as  the 
ballistic  damage  requirements.  The  design  weight  for  the  metal 
wingbox  increased  more  than  the  composite  design  for  the  following 
reason.  In  sizing  the  metal  wingbox,  a  comparison  was  made  to  the 
baseline  A-7D  outer  wing.  The  skins  of  the  baseline  wingbox  were 
thicker  in  some  regions  than  the  ADDRESS  output  and  the  Vought 
input  (from  a  NASTRAN  model) .  The  metal  skins  were  sized  according 
to  this  information  to  assure  that  the  wingbox  met  all  static,  fatigue 
and  fracture  requirements,  A  comparison  of  the  computer  runs  for  the 
metal  wingbox  shows  that  the  majority  of  the  increase  from  130.7  lbs. 
to  160.55  lbs.  was  due  to  the  increase  in  skin  weights. 


* 
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SECTION  7 

COMMENTS  AND  CONCLUSIONS 

An  analytical  procedure  has  been  developed  and  demonstrated 
which  provides  optimal  design  information  for  lifting  surfaces 
subjected  to  strength ,  stiffness ,  and  damage  requirements .  The 
optimality  criterion  approach  provides  a  framework  in  which  effi¬ 
cient  reanalysis  techniques  can  be  used  to  handle  the  multiple 
requirements.  An  examination  of  the  optimal  design  which  includes 
the  effect  of  damage  shows  that  the  ADDRESS  program  provides  alter¬ 
nate  load  paths  around  damaged  structural  members. 

Execution  times  have  been  reduced  by  33%  to  50%  compared  to 
similar  computer  runs  using  0PSTATC0MP  or  DANALYZ .  The  ADDRESS 
program  because  of  its  modularity  provides  a  framework  within 
which  new  program  features  can  be  easily  added. 

Computer  output  has  been  used  by  an  airframe  company  (Vought 
Corporation)  to  establish  viable  designs  for  given  damage  con¬ 
ditions.  Even  with  a  fairly  coarse  finite  element  model,  useful 
information  can  be  obtained  which  assists  the  designer  in  arriving 
at  lower  weight,  damage  tolerant  designs. 

One  disadvantage  of  the  method  appears  to  be  that  resulting 
designs  obtained  from  ADDRESS  have  locally  increased  member  sizes 
near  the  damage  zone.  Some  type  of  redistribution  technique  is 
needed  to  compensate  for  this  problem.  One  approach  is  to  run  the 
ADDRESS  program  for  a  wide  variety  of  damage  conditions  and  then 
take  as  the  candidate  design  the  envelope  of  the  maximum  member 
sizes  so  obtained.  Another  approach  would  be  to  impose  a  low  level 
damage  to  every  member  and  then  optimize  for  this  single  uniform 
damage  condition. 

In  general  a  more  organized  approach  is  needed  for  the  selec¬ 
tion  of  the  most  critical  damage  conditions.  The  capability  is 
available  within  ADDRESS  to  treat  this  problem  in  a  limited  manner. 
If  is  possible  to  damage  any  individual  element  and  determine  by 
reanalysis  the  response  of  the  structure.  The  most  critical 
elements  or  groupings  of  elements  could  thus  be  identified  and 
appropriate  optimization  runs  then  performed.  The  ADDRESS  program 
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is  currently  limited  in  the  number  of  damage  cases  it  can  handle 
and  would  need  some  modifications  to  be  used  in  the  above  manner. 

The  problem  also  exists  that  it  is  difficult  to  convert  a 
given  ballistic  damage  condition  into  a  single  mass  or  stiffness 
reduction  factor  for  the  ADDRESS  program.  The  bar,  shear  panel, 
and  membrane  structural  model  has  limitations  in  representing  com¬ 
plicated  ballistic  damage  effects.  On  the  other  hand  it  is 
j,  probable  that  far  from  the  local  damage,  the  member  sizing  is 

independent  of  the  details  of  the  damage.  Further  studies  on  a 
*  component  level  would  be  needed  to  establish  a  sound  basis  for  the 

selection  of  reduction  factors  to  be  applied  to  the  bar,  shear 
panel,  and  membrane  model.  Mathematical  programming  techniques 
applied  at  the  component  level  in  a  separate  analysis  and  optimi¬ 
zation  program  could  provide  useful  input  the  the  less  sophi¬ 
sticated  overall  model  which  is  resized  by  optimality  criteria. 

Further  work  is  also  needed  to  quantify  the  effects  of  bal¬ 
listic  damage  on  the  applied  loads.  From  an  aerodynamic  point  of 
view,  this  load  redefinition  could  be  quite  complex  due  to  non- 
uniform  airflow  through  or  near  the  ballistic  damage. 

Other  constraints  need  to  be  considered  in  the  damage  optimi¬ 
zation  problem.  A  constraint  on  the  ratio  of  bending  to  torsional 
frequencies  for  a  high  aspect  ratio  wing  could  be  quite  important 
from  an  aeroelastic  stability  standpoint.  Constraints  on  diver¬ 
gence  and  flutter  speeds  need  to  be  considered  together  with  the 
strength  requirements. 
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